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CHAPTER  1 
Introduction 

1.1  Objective 

The  objective  of  this  study  is  to  investigate  the  propagation  of  acoustic  waves 
in  contiguous  concentric  cylinders  of  solid  and  fluid  media.  Wavenumber-frequency 
and/or  transform  methods  are  applied  to  several  configurations  of  interest  to  de¬ 
termine  the  characteristics  of  acoustic  wave  propagation  as  functions  of  the  spatial 
variables  and  the  properties  of  the  media.  In  the  experimental  work,  joint  time- 
frequency  distribution  functions  such  as  the  spectogram  are  used  to  study  the 
propagation  of  transient  disturbances  in  dispersive  systems.  The  study  is  limited 
to  systems  of  infinite  extent  in  the  axial  direction. 

1.2  Background 

There  are  many  practical  situations  in  which  it  is  desireable  to  understand  the 
characteristics  of  the  propagation  of  acoustical  energy  in  physical  systems  which 
may  be  modeled  in  a  cylindrical  coordinate  system.  For  example,  in  noise  control 
problems,  the  objective  is  to  reduce  the  transmission  of  sound  from  one  region  into 
another.  Cases  of  interest  include  aircraft  fuselages  and  submarine  hulls  in  which 
turbulent  boundary  layer  noise  is  generated  by  the  motion  of  the  platform  through 
an  external  fluid  medium.  Another  case  is  noise  generated  by  fluid  flow  in  pipes 
and  ducts.  Another  problem  of  interest  relates  to  sonar  transducers,  arrays,  and 
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domes  of  cylindrical  geometry,  in  which  it  is  desired  to  enhance  the  transmission 
of  propagating  acoustic  waves  in  the  external  medium  while  discriminating  against 
evanescent  waves  generated  by  turbulence  or  structural  vibration. 

There  is  an  extensive  body  of  literature,  including  both  textbooks  and  mono¬ 
graphs  in  technical  journals,  that  deals  with  the  problems  considered  in  this  dis¬ 
sertation.  Specific  subjects  of  interest  include  the  propagation  of  acoustic  waves  in 
ducts,  acoustic  radiation  in  an  infinite  fluid  medium,  the  transmission  of  acoustical 
waves  in  layered  media,  and  the  vibration  of  shells. 

In  problems  which  involve  fluid-structure  interaction,  the  most  complex  task  is 
analysis  of  the  motion  of  the  structure.  Thus  it  is  customary  to  first  consider  the 
"vibration  of  the  structure  in  vacuo  and  then  to  attempt  to  determine  the  effects  of 
coupling  to  the  fluid. 


1.3  Organization  of  the  Study 

The  general  schema  of  the  models  studied  in  this  work  is  outlined  in  Figure 
1-1.  The  general  category  of  cylindrical  models  branches  into  fluid  media  and 
solid  media.  Under  fluid  media,  models  are  developed  for  the  fluid  cylinder,  a 
fluid  tube,  and  the  external  fluid  medium.  The  branch  comprising  elastic  media  is 
subdivided  to  include  models  which  may  be  represented  by  thin  shell  theory,  and 
those  for  which  elasticity  theory  is  required.  This  latter  category  includes  thick 
shells,  concentric  elastic  shells,  and  the  external  solid  medium.  Composite  system 
models  are  constructed  by  combining  the  elementary  models  of  system  components 
through  the  transfer  matrix  approach. 
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In  this  introductory  chapter,  we  describe  the  transfer  matrix  method  which 
was  suggested  by  previous  work  on  the  propagation  of  plane  waves  through  layered 
media.  This  is  followed  by  a  discussion  of  wavevectors  in  cylindrical  coordinates, 
and  the  synthesis  of  cylindrical  waves  from  plane  waves. 

Chapter  2  begins  with  a  discussion  of  the  analytical  approach  followed  for 
the  study  of  fluid  media.  These  methods  are  applied  to  the  fluid  cylinder,  the 
fluid  tube,  and  the  external  medium.  For  each  of  these  components,  the  radial 
spectral  impedance  is  computed.  The  dispersion  equation  is  developed  and  used 
to  study  the  propagation  of  free  waves.  The  transfer  matrix  for  each  component 
is  determined. 

In  Chapter  3  we  develop  the  transfer  matrix  for  a  thin  elastic  cylinder,  using 
Donnell’s  formulation  (Donnell  1933;  Krauss  1967;  Junger  and  Feit  1986)  Disper¬ 
sion  of  harmonic  waves  for  the  axysimmetric  mode  is  studied  in  the  frequency 
domain.  An  analysis  of  low  frequency  axisymmetric  longitudinal  modes  of  a  har¬ 
monically  driven  finite  length  cylinder  is  presented.  The  transient  response  of  a 
semi-infinite,  radially  excited  elastic  cylinder  is  studied  using  Laplace  and  Fourier 
transforms  in  the  wavenumber  and  frequency  domains. 

The  models  developed  in  the  first  three  chapters  are  combined  in  Chapter  4  to 
describe  several  composite  systems  of  interest. 

Chapter  5  describes  the  results  of  transient  measurements  on  a  thin  brass 
cylinder.  These  results  are  compared  to  the  predictions  made  in  Chapter  3. 
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1.4  Transfer  Matrix  Approach 

The  transfer  matrix  approach  followed  in  this  study  is  similar  to  that  employed 
by  pervious  authors  in  analyzing  the  propagation  of  plane  waves  in  layered  media. 
Thomson  (1950)  introduced  the  matrix  method  to  study  the  transmission  of  plane 
elastic  waves  through  any  number  of  parallel  plates  of  solid  or  fluid  material.  Shaw 
and  Bugl  (1969)  employed  matrix  methods  to  study  harmonic  plane  waves  in 
viscoelastic  media,  using  complex  and  frequency  dependent  Lame  constants  to 
describe  the  viscoelastic  behavior.  Folds  and  Loggins  (1977)  extended  this  work 
to  include  calculations  of  the  transmission  and  reflection  coefficients.  They  report 
experimental  work  which  agrees  well  with  the  theory.  Stepanishen  and  Strozeski 
(1982)  employed  transfer  matrix  and  Fourier  transform  methods  to  analyze  the 
response  of  layered  viscoelastic  media  to  a  time  dependent  plane  wave  excitation. 

In  this  study  models  are  constructed  for  each  component  of  the  concentric 
cylindrical  system.  Typical  components  include,  for  example,  a  fluid  cylinder,  a 
fluid  tube,  an  elastic  shell,  and  the  external  medium. 

For  each  component,  the  governing  equations  for  the  field  variables  in  terms  of 
time  and  spatial  coordinates  are  converted  via  a  two  dimensional  Fourier  transform 
into  a  representation  in  wavevector-frequency  space.  The  transfer  function  of  the 
component  is  then  formulated  in  terms  of  the  state  variables  at  the  radial  interfaces. 
The  resulting  matrix  model  is  two  by  two  (two  port)  for  fluid-solid  and  fluid-fluid 
interfaces,  for  which  the  boundary  conditions  are  continuity  of  pressure  and  radial 
velocity.  For  contiguous  elastic  solids  in  which  shear  stresses  exist,  continuity  of  the 
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three  components  of  displacement  and  stress  results  in  a  six  port  transfer  matrix. 

A  system  model  for  any  combination  of  components  is  constructed  by  cas¬ 
cading  the  component  transfer  matrices  to  obtain  an  overall  transfer  function. 
Specification  of  an  appropriate  spectral  excitation  function  in  one  region  of  the 
system  model  enables  the  computation  of  state  and  field  variables  throughout  the 
system.  Non-axisymmetric  excitations  may  be  represented  using  a  Fourier  series 
expansion  in  terms  of  a  discrete  circumferential  wavenumber. 

Thus  the  general  problem  of  interest  comprises  a  number  of  concentric  circular 
cylinders  imbedded  in  an  infinite  fluid  medium.  This  system  is  excited  by  a  pressure 
or  velocity  function  specified  at  one  of  the  radial  interfaces.  Figure  1-2  indicates 
the  coordinates  of  a  general  cylindrical  system  with  m  components  in  which  the 
center  component  is  a  cylinder  of  radius  Ui,  surrounded  by  a  tube  of  inner  radius 
ai  and  outer  radius  a2,  etc.  terminating  in  the  (m  —  1  )th  tube  of  inner  radius  dj 
and  outer  radius  am,  which  marks  the  boundary  with  the  external  medium.  The 
general  cylindrical  coordinates  are  r,  0,  and  2,  representing  the  radius,  azimuthal 
angle,  and  cylindrical  axis,  respectively. 

The  two-port  model  for  a  fluid  tube  is  shown  in  Figure  1-3.  If  this  is  the  j  th 
component  in  the  system,  with  inner  radius  and  outer  radius  a3  ,  the  state 

variables  of  pressure  and  velocity  V(j_i)  at  the  inner  radius  are  related  to 

the  pressure  Pj  and  velocity  V3  at  the  outer  radius  by  a  transfer  matrix  [Tj]  ,  viz. 


- 1 

_  J 

pu- 1) 

=  Pi] x 

1 

S’ 

_ 1 

1 

1 

f— * 

1 _ 

The  outermost  cylinder  is  surrounded  by  a  fluid  presumed  to  extend  to  infinity 
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in  the  radial  and  axial  directions.  Application  of  the  radiation  condition  that  the 
pressure  must  decrease  to  zero  at  infinity  leads  to  equations  for  the  specific  acoustic 
impedance  at  the  last  radial  interface  at  r  =  am. 

A  further  requirement  is  that  all  concentric  layers  must  exhibit  the  same  pres¬ 
sure  dependence  in  the  axial  direction  as  the  excitation. 

For  example,  if  the  excitation  is  a  traveling  wave  represented  by 

pn(z,  6,  a,j )  =  Pe  cos (n9)  exp  (i  ( kzz  —  cut)) 

where  kz  is  the  axial  wavenumber  of  the  excitation,  u  is  the  angular  frequency, 
t  represents  time,  and  n  is  the  circumferential  mode  number,  then  each  component 
of  the  system  must  have  the  same  variation  with  kz  and  n.  This  condition  results 
from  the  continuity  of  normal  particle  velocity  and  pressure  across  the  interface 
between  contiguous  layers. 

1.5  Wavevectors 

The  concept  of  a  wavevector,  also  known  as  a  vector  wavenumber  (Elmore  & 
Heald,  1985)  is  highly  useful  in  describing  complex  wave  phenomena.  We  define 
the  wavevector  k  for  a  monochromatic  wave  as  a  quantity  having  the  magnitude 
of  a  wavenumber,  |k|  =  j  and  the  direction  of  a  unit  vector  which  is  normal  to 
the  wavefront,  or  surfaces  of  constant  phase.  For  example,  a  propagating  wave  of 
frequency  uj  may  be  expressed  compactly  as 

p(r,t)  =  Re{P0e<rk-^}  (1.1) 

where  PQ  represents  the  complex  wave  amplitude  and  r  is  the  spatial  vector.  In 
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rectangular  coordinates,  r  =  exx+eyy+ezz  where  ex,ey,ez  are  unit  vectors  along 
the  x ,  y,  and  z  axes,  respectively.  In  this  formulation,  surfaces  of  constant  phase  are 
described  by  r  •  k  =constant.  In  a  region  with  no  acoustic  sources,  equation  (1.1) 
must  satisfy  the  reduced  (Helmholtz)  wave  equation,  ^V2  +  |k|2)  P  =  0;  application 
of  this  to  (1.1)  yields  the  result  |k|2  =  =  A:2  +  A;2  +  k\  .  The  domain  of 

the  wavevector  k  is  thus  described  by  a  sphere  of  radius  |k|  ,  as  indicated  in 
Figure  1-4.  Propagating  acoustic  waves  are  represented  by  points  on  the  surface 
of  this  sphere.  Wavevectors  falling  within  the  surface  of  the  sphere,  where  |k|  < 

.  describe  supersonic  waves.  Wavevectors  exterior  to  the  surface  of  the  sphere 
describe  evanescent  sub-sonic  waves.  Note  that  the  quantities  and  ^  are 

the  direction  cosines  of  k  with  respect  to  the  x,y.  and  z  axes.  We  identify  cx  = 
r- .  cu  =  t-  ,  and  cz  =  as  the  trace  velocities  of  the  wave  along  the  respective 

K-x  *  fZy  Kz 

axes.  Similarly,  since  k  =  ^  where  A  is  the  wavelength  of  the  propagating  wave, 
the  trace  wavelengths  are  expressed  as  A^  =  — ■,  Xy  —  |^,  and  A 2  =  f~* 

We  note  that  as  the  wavevector  k  approaches  perpendicularity  with  respect 
to  one  of  the  coordinate  axes,  the  trace  velocity  and  wavelength  tend  to  infinity 
along  that  axis.  This  concept  is  useful  in  understanding  the  cutoff  phenonmenon 
which  occurs  under  certain  conditions  in  a  propagating  wave  which  is  constrained 
by  some  physical  boundary,  such  as  within  a  waveguide. 

If  p  (r,  t )  describes  a  complex  wave  field  for  which  the  principle  of  superposition 
applies,  we  may  define  a  wavevector  spectrum  in  terms  of  a  three  dimensional 
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spatial  Fourier  transform, 

OO 

P  (k,  *)  =  JJJ P  (r,  *)  e~ik  T 

—  OO 

in  which  p(k, t)  =  p{kx,ky,kz,t)  represents  the  spectral  components  of  the 
wavefield  in  the  wavevector  domain  in  the  same  manner  that  the  Fourier  transform 
of  a  time  dependent  function  depicts  the  spectral  composition  of  the  function  in 
the  frequency  domain. 

In  this  study  we  are  primarily  concerned  with  wave  propagation  in  cylindrical 
configurations  with  coordinates  ( r,8,z )  and  unit  vectors  (er,  eg,  ez).  In  this  no¬ 
tation,  the  wavevector  may  be  written  kcyi  =  erkr  +  egkg  +  ezkz.  A  vector  krect  in 
rectangular  coordinates  may  be  transformed  into  a  cylindrical  wavevector  k cyi  in 
cylindrical  coordinates  by  means  of  an  orthogonal  transformation  matrix  [5]  which 
has  the  following  properties  (Hughes  and  Gaylord  1964): 

1.  The  inverse  of  [S']  is  equal  to  the  transpose  of  [S']  ,  [S]_1  =  [S]T ,  and 
[S]  ■  [S]  1  =  [I],  the  identity  matrix. 

2.  The  determinant  of  [S]  =  ±1. 

The  transformation  equations  are:  k^/  =  [S]  •  k rect,  and  kreci  =  [S]T  •  k cy[. 
where 

cos  6  sin  9  0 
[PI  =  —  sin  8  cos  8  0 


0 


0  1 
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and 

cos  9  —  sin  9  0 
[£]T  —  sin  9  cos  9  0 

0  0  1 

Applying  this  to  the  present  case,  the  components  of  the  cylindrical  wavevector 
in  terms  of  the  elements  of  a  rectangular  wavevector  are: 

kT  =  kx  cos  9  +  ky  sin  9;  ke  =  —  kx  sin  9  +  ky  cos  9;  and  kz  =  kz. 

The  inverse  transformation  results  in: 

kx  =  kT  cos  9  —  ke  sin  9\  ky  =  kT  sin  9  +  ko  cos  9 ;  and  kz  —  kz. 

Propagating  cylindrical  waves  must  satisfy  the  equation 

|k|2  =  kl  +  k%  +  fcf  (1-2) 

Now  let  ke  =  0.  For  this  case,  Eq.  (1.2)  defines  a  cone  in  k  -  u  space, 
Figure  1-5.  Points  on  the  surface  of  this  cone  correspond  to  propagating  acoustic 
waves.  A  plane  at  some  particular  frequency  ujq  which  is  parallel  to  the  u  =  0 
plane  intersects  the  cone  and  defines  a  circle  of  radius  ko.  Figure  1-6.  The  two- 
dimensional  wavevector  k  =  erkr  +  ezkz  represents  propagating  waves  when  |k|  = 
A:0 .  supersonic  waves  when  k  is  inside  the  circle,  and  sub-sonic  evanescent  waves 
when  k  is  outside  the  circle. 


A  construction  which  is  useful  in  the  description  of  cylindrical  waves  comprises 
the  synthesis  of  such  waves  from  an  infinite  number  of  plane  waves.  (Lamont  1950; 
Redwood  1960) 
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Consider  a  plane  wave  of  wavenumber  k  propagating  in  a  rectangular  ( x ,  y,  z) 
coordinate  system,  with  wavevector  k  making  an  angle  (p  with  respect  to  the  z 
axis,  and  oriented  so  that  there  is  no  component  in  the  y  direction,  i.  e.,  ky  =  0. 
Such  a  wave  may  be  described  by 

p  (x,  y,  z,  t )  =  Re  {i4e<(Jfe-!B+**z"wt)} 
in  which  kx  =  k  sirup  and  kz  =  k  cos  <p. 

Now  we  hold  the  angle  (p  constant  and  sum  up  an  infinite  number  of  such  plane 
waves  by  rotating  the  x,y  axes  around  the  z  axis.  We  consider  a  point  at  (r,  6,0) 
in  the  plane  z  =  0,  Figure  1.7.  The  x  axis  is  shown  an  an  angle  ip  with  respect  to 
the  polar  coordinates  (r,  0,0)  and  thus  x  —  r  cos  ip.  Now  let  the  amplitude  of  the 
elementary  wave  associated  with  the  point  (r,  6,  z)  be  represented  by 

dp  =  Re  {Acosn  ( 6  +  ip)  e^ucosi>+kzz-u)  ^  j 

The  total  pressure  at  (r,  6,  z)  is  obtained  by  integrating  ip  over  the  range  0  < 
ip  <  27 r, 

U27T 

A  COS  n  {9  +  Ip)  e^rsmtcos^jikzcost-wt)  ^ 

Using  the  identity  (Gradshteyn  and  Ryzhik  1980) 

1 ~n  P2  7T 

Jn  (u)  =  —  gi(u cosip)  CQUn^  dxp 

2sk  Jq 

and  substituting  u  =  kr  sin  0,  we  obtain 


V  (r,  e,  z,t)  =  Re  {27rinAcos  n6  Jn  (kr  sin  <j>)  J  (1.3) 
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Making  the  substitutions  kT  =  k  sin0  and  kz  =  k  cos  0  in  Eq.  (1.3)  yields 

p  (r,  6,  z,t)  =  Re  |27rin^4cos  nO  Jn  (krr)  e*(^2_a;b|  (1-4) 

which  is  recognized  as  one  solution  of  the  cylindrical  wave  equation,  and  kr  and 
kz  are  components  of  the  cylindrical  wave  vector.  For  the  axisymmetric  (n  =  0) 
case,  Eq  (1.4)  specializes  to 

p  (r,  z,  t )  =  Re  {2t rA  J0  ( kTr ) 

This  describes  a  wave  diverging  from  the  origin  in  a  cone  of  half-angle  (p ,  as 
indicated  in  Figure  1.8.  If  the  amplitude  constant  A  is  real,  this  may  be  written 


p  ( r ,  z,  t )  =  27T.A  J0  (krr)  cos  ( kzz  -  cut) 


Radial:  vT  (r,  z,  t )  =  —  27t^  sin  <p  Ji  ( kTr )  sin  ( kzz  —  cut) 
and  Axial:  vz  (r,  z,  t)  =  27r^  cos  0  Jo  (krr)  cos  ( kzz  —  wt) 

Since  the  axial  velocity  and  the  pressure  are  in  phase,  energy  is  propagating 


in  the  axial  direction.  The  radial  velocity  is  orthogonal  to  the  the  pressure. 


12 


CHAPTER  2 
Fluid  Media 


2.1  Analytical  Approach 


In  this  chapter  we  develop  the  transfer  matrices  for  three  cylindrical  fluid 
components  of  interest:  a  fluid  cylinder,  a  fluid  tube,  and  the  external  fluid.  Wave 
propagation  in  these  components  is  studied  for  limiting  boundary  conditions  (rigid 
and  soft)  as  a  prelude  to  understanding  the  more  complicated  situation  which 
obtains  when  composite  (cascaded)  systems  are  considered. 

For  all  fluid  models  we  shall  assume  the  fluid  is  inviscid  and  that  the  amplitudes 
of  the  field  variables  (pressure  and  particle  velocity)  are  sufficiently  small  that  the 
linearized  wave  equation  is  applicable.  Thus  we  are  concerned  with  solutions  to 


Vzp 


■2_  1  d2P 


c2  dt2 


where  p  is  the  acoustic  pressure  in  the  fluid,  c  is  the  acoustic  wave  speed  in  the 
fluid,  and  V2  is  the  three-dimensional  Laplacian  operator.  The  particle  velocity  v 
is  related  to  the  pressure  via  the  linearized  momentum  equation, 


Vp  =  -p0 


dv 

dt 


where  po  is  the  equilibrium  (static)  density  of  the  fluid. 


2.1.1  Wave  Equation  in  Cylindrical  Coordinates 

Consider  the  cylindrical  coordinate  system  shown  in  Figure  2.1,  in  which  r,z  and 


6  are  the  radial,  axial,  and  angular  variables,  respectively;  and  er,  e$,  ez  are  the 
corresponding  unit  vectors.  The  wave  equation  for  the  pressure  is 


Id/  dp\  1  d2p  d2p  _  1  d2p 
r  dr  1 T  dr)  r2  dd2  dz 2  c2  dt 2 


(2.1) 


and  the  momentum  equation  is 


dp  1  dp  dp 
e'fr+eerd§+e‘Tz  = 


‘Po- 


dv 

dt 


(2.2) 


Two  approaches  in  solving  the  cylindrical  wave  equation  (2.1)  are  the  separa¬ 
tion  of  variables  and  the  use  of  integral  transforms ,  such  as  the  Fourier  or  Laplace 
transforms. 

In  the  separation  of  variables  method,  the  solution  is  assumed  to  consist  of  the 
product  of  independent  functions  of  r,  z,  6  and  t,  such  as 


p(r,9,z,t)  =  R(r)e(e)Z(z)T(t)  (2.3) 

If  harmonic  time  dependence  of  the  form  T(t )  =  exp (—iut)  is  assumed,  the 
reduced  (harmonic)  Helmholtz  equation  is  obtained, 

(V2  +  k2)p(r,e,z)=  0  (2.4) 


where  is  the  angular  frequency  and  k  =  —  is  the  acoustic  wavenumber. 

Operating  on  Eq.2.3  using  Eq.2.4,  and  suppressing  the  time  dependence,  three 
independent  ordinary  differential  equations  are  obtained: 


In  these  equations,  k2 .  kr  and  n  are  separation  constants,  related  by  the  ex¬ 
pression  k2  —  k2  -I-  k2  +  j  . 

Solutions  to  2.5  may  be  either  exponential,  Z  (z)  =  Ae(±ikzZ\  representing 
traveling  waves,  or  sinusoidal,  ,  for  standing  waves. 

Since  0  ( 6 )  repeats  periodically  with  period  2tt  ,  that  is,  ©  (0)  =  ©  {0  +  27t), 
n  must  be  an  integer,  and  solutions  of  2.6  must  be  combinations  of  trigonometric 
functions,  such  as  ©  ( 9 )  =  B\  cos  (nO)  +  £>2  sin  (n8) . 

Equation  2.7  is  recognized  as  Bessel’s  equation  of  integer  order  n,  solutions  to 
which  comprise  various  linear  combinations  of  Bessel  functions  of  argument  krr. 
such  as: 

Rn  {r)  —  Jn  ( kTr )  +  iC2Yn  ( krr ) 

or  alternatively, 

Rn  (r)  =  DXH^  ( Ayr )  +  D2H. <2)  (krr) 

In  these  equations,  kz  may  be  identified  as  an  axial  wavenumber,  and  kr  as  a 
radial  wavenumber.  At  a  particular  value  of  the  radial  coordinate,  r  =  a,  we  may 
define  kn  =  **,  a  discrete  circumferential  wavenumber.  By  use  of  the  momentum 
equation  at  the  radial  interfaces  r  =  <2^0,2,  ...arn  for  a  particular  cylindrical  fluid 
component,  and  recognizing  that  kz  and  n  must  be  identical  for  each  component 
in  the  system  to  insure  the  continuity  of  the  radial  particle  velocity  across  each 
radial  interface,  the  various  coefficients  A,  B,C  etc  may  be  evaluated. 
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2.1.2  Integral  Transform  Approach 

A  major  advantage  of  integral  transforms,  such  as  the  Fourier  or  Laplace  trans¬ 
forms,  is  the  reduction  of  differential  equations  to  algebraic  equations.  In  the  solu¬ 
tion  of  transient  problems,  the  Laplace  transform  is  advantageous  since  it  permits 
the  initial  conditions  to  be  included  in  the  formulation  of  the  problem. 

To  illustrate  the  application  of  integral  transforms  we  will  construct  a  Fourier 
transformed  version  of  the  cylindrical  wave  equation,  transforming  on  the  time  t, 
and  the  spatial  coordinates  z  and  6. 

For  the  continuous  variables  z  and  t  we  define  the  forward  transform 

oo 

p(kz,oj)  =  JJp(z,t)e-^t+k‘z)  dt  dz  (2.8) 

—  OO 

The  inverse  transform  is  defined  as 


i  00 

pCM)  =  j^Jf  iu, 


dk7 


Since  the  circumferential  variation  is  periodic,  we  may  expand  p  (0)  in  a  Fourier 


series,  viz. 


p(e)  =  S 


p  (n)  e 


ind 


27 r 


(2.9) 


Multiplying  both  sides  of  Eq.  (2.9)  by  e~inU>  and  integrating  over  the  range 
— 7T,  7r  we  have 

r  p(e)e-im6dd=  T  r  i^-e^n-m)ede  (2.10) 

J-*  n=— oo  J-*  2?r 

From  the  orthogonality  of  trigonometric  functions,  the  integral  on  the  right 
hand  side  of  Eq.  (2.10)  is  zero  for  all  values  of  n  ^  m  and  therefore 


p(n)=  r  p(e)e-indd9 

J  -  7T 


(2.11) 
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Equations  (2.9)  and  (2.11)  define  a  transform  pair  on  the  circumferential  co¬ 
ordinate  6  . 


Transforms  of  Derivatives  Differentiating  Eq.  (2.9)  m  times  yields 


_  (jB)-  f  me 

dQm  v  }  ^  2tt 

n— — oo 


ind 


and  transforming  both  sides  of  this  equation  gives  for  the  transform  of  the  deriva¬ 


tive: 


T 


dmp  (8) 


dQ« 


=  (' in)mp(d ) 


The  derivatives  of  the  transforms  of  the  continuous  variables  2  and  t  are  ob¬ 


tained  by  differentiating  Eq.  (2.8),  viz 

'dmp(tY 


T 


dtT 


=  (-iu)mp(u>) 


and 


T 


dmp  ( z ) 
dzm 


(~ikz)mp(kz) 


2.1.3  Transformed  Cylindrical  Wave  Equation 

Applying  these  transforms  to  the  cylindrical  wave  equation  (2.1)  results  in  an 
ordinary  differential  equation  in  r  and  in  p(kz,u,n) , 


d?p  l  dp  (  2  2  n2\  _ 

d?  +  rTr  +  \k 


(2.12) 


which  is  identical  in  form  to  Eq.  (2.7)  obtained  by  the  method  of  separation  of 
variables.  Thus  solutions  of  (2.12)  will  consist  of  combinations  of  Bessel  functions 
of  argument  krr  and  coefficients  determined  by  the  boundary  conditions  which 


apply  to  the  model  of  interest. 
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2.2  Fluid  Cylinder 


2.2.1  Introduction 

The  propagation  of  acoustic  waves  in  a  confined  fluid  has  received  extensive  study 
in  the  literature.  Redwood  (1906)  discusses  continuous  and  transient  (pulsed) 
propagation  in  rectangular  and  cylindrical  fluid  waveguides  with  both  free  and  rigid 
boundaries.  Morse  and  Ingard  (1968)  devote  a  chapter  to  ”  Sound  Waves  in  Ducts 
and  Rooms”  in  which  both  plane  wave  and  higher  order  modes  are  considered. 
Their  analysis  includes  an  intermediate  boundary  condition  in  which  the  stiffness 
and  mass  of  the  duct  walls  are  considered.  A  similar  analysis  is  presented  by 
Junger  and  Feit  (1986). 

Anderson  and  Barnes  (1953)  studied  the  dispersion  of  a  damped  sinusoidal 
ultrasonic  pulse  in  a  water-filled  thin  walled  cylinder  lined  with  rubber  to  create 
a  free  radial  boundary.  Their  experimental  results  agreed  with  the  theoretical 
dispersion  curve  for  the  first  axisymmetric  (0,1)  mode.  They  also  developed  an 
expression  for  the  transient  response  as  a  function  of  time  and  axial  distance  for  an 
impulsive  (delta  function)  velocity  input  and  related  this  result  to  the  dispersion 
predicted  from  group  velocity  concepts. 

The  dispersion  of  a  rectangular  narrow  band  pulse  (tone  burst)  in  a  rectan¬ 
gular  fluid  waveguide  with  free  boundaries  was  studied  by  Proud,  Tamarkin,  and 
Kornhauser  (1956).  They  measured  the  phase  velocity  as  a  function  of  frequency 
for  the  (1,1)  and  (2,1)  modes  of  the  waveguide.  The  pulse  envelopes  were  measured 
at  various  distances  along  the  waveguide  and  compared  with  the  predictions  from 


18 


a  theoretical  model  by  Pearson  (1953). 

Jacobi  (1949)  made  an  extensive  theoretical  and  experimental  study  of  the 
propagation  of  acoustic  waves  in  fluid  cylinders.  His  models  included  a  fluid  cylin¬ 
der  with  rigid  and  free  walls,  a  fluid  cylinder  imbedded  in  an  infinite  fluid,  a  fluid 
cylinder  with  fluid  walls,  and  a  fluid  cylinder  within  a  thin  elastic  cylinder.  Exper¬ 
imental  results  were  obtained  for  the  fluid  cylinder  with  free  boundaries  and  for 
the  fluid  cylinder  in  a  thin  shell. 


2.2.2  Physical  Description 

Consider  a  cylinder  of  fluid  of  radius  a,  sound  speed  c,  and  density  p  extending 
infinitely  in  the  axial  (z)  direction.  The  solution  of  the  transformed  wave  equation 


(2.12)  appropriate  for  this  case  is 


p(n,kz,u,r)  =  p  (n,  k2 , u)  Jn  ( kTr ) , 

where  k%  =  —  k\  —  .  (The  Neumann  function  Nn  ( kTr )  is  excluded 

because  of  the  singularity  at  r  =  0.) 

The  transformed  momentum  equation  is 

dp  (in)  _  ^  ,  _  .  ,  _ 

er—  +  e# - p  +  ez  (xkz)  p  =  p  (zu>)  v  =  ipckv 

Application  of  the  momentum  equation  to  the  transformed  wave  equation  re¬ 
sults  in  the  three  components  of  particle  velocity: 

Axial: 

kz 
pck 


p  Jn  ( krr ) 


V 
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Circumferential: 


Vo  = 


n 


pckr 


p  Jn  (hr) 


Radial: 


(2.13) 


where  the  prime  denotes  the  derivative  with  respect  to  the  argument  (kTr ) . 
These  relationships  will  be  used  to  develop  the  radial  spectral  impedance  of 
the  fluid  cylinder,  and  to  study  the  propagation  of  free  waves  in  a  fluid  cylinder. 

2.2.3  Radial  Spectral  Impedance 

We  define  a  radial  spectral  impedance  as  the  ratio  of  the  transformed  pressure  to 
the  transformed  radial  component  of  velocity; 

^  /  ,  x  P  ipck  Jn  (hr) 

z(^M  =  Wr=  j,(M- 

We  now  specialize  this  to  the  axisymmetric  (n  —  0)  mode.  The  justification  for 
this  simplification  is  that  any  arbitrary  circumferential  excitation  function  may  be 
expanded  in  a  Fourier  series  as  indicated  in  Eq.  (2.9).  Solutions  would  then  consist 
of  a  summation  of  terms  representing  the  harmonic  components  of  the  excitation. 

For  the  axisymmetric  case,  J'Q  ( krr )  =  -kr  Jj  (krr)  and  the  spectral  impedance 
of  the  fluid  cylinder  is 

z(u  k  .-ipckMkrr) 

A;r  J\  (krr) 

At  the  surface  of  the  cylinder,  r  =  a,  we  define  the  state  variables  of  pressure 
Pa  and  Va  in  terms  of  the  spectral  impedance,  viz: 

Pa  —ipck  Jo  (ha) 


Za 


Va  h  J\  (ha) 


(2,14) 
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The  notation  may  be  simplified  somewhat  by  using  the  definition  for  the  mod¬ 
ified  quotients  of  Bessel  functions  (Onoe  1958).  Therefore  we  define  Jn  (kra)  as 

(kT<i}  Jn+l  ( kro ;) 

Jn  ( krCL 

This  function  is  plotted  in  Figure  2.2  for  n  =  0  and  values  of  the  argument 
from  0  to  15. 

With  this  definition,  the  radial  spectral  impedance  of  Eq.  (2.14)  is  written  as 

—ipcka 

\Jo  ( kra ) 

To  study  the  spectral  impedance  of  the  fluid  cylinder,  consider  an  excitation 
consisting  of  a  spectrally  pure  traveling  wave  of  pressure  of  arbitrary  frequency  u> 
and  axial  wavenumber  kz  applied  at  the  boundary  r  =  a.  This  condition  could  be 
physically  realized,  for  example,  by  a  membrane  or  thin  elastic  cylinder  enclosing 
the  fluid  and  driven  by  some  external  force. 

The  waves  in  the  fluid  cylinder  must  satisfy  the  dispersion  equation  k%  = 
(y  j  —  k^.  Thus  as  lu  and  kz  take  on  all  possible  values  the  radial  wavenumber  kr 
may  be  real,  zero,  or  imaginary. 

Case  1  -  kr  Real,  Propagating  Axial  Wave  For  k  >  kz,  kr  is  real,  and  kr  = 
yjk2  —  fc|.  The  normalized  spectral  impedance,  Za/pc,  is  plotted  in  Figure  2.3  as 
a  function  of  the  argument  kra  for  ka  =  1. 

The  spectral  impedance  is  purely  reactive,  alternating  between  stiffness  and 
mass  type  reactances  as  the  radial  wavenumber  kra  increases.  With  kz  real,  waves 
may  propagate  along  the  z  axis  of  the  fluid  cylinder. 
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Case  2  -  No  Axial  Wave  Propagation  Now  if  we  let  kz  =  0,  (i.  e.  there  is  no  axial 
propagation  in  the  fluid  cylinder)  then  kr  =  k,  and 


— ipc  J0  ( ka ) 
Jj  (ka) 


(2.15) 


For  low  frequencies  (ka  <Cl,  or  u  c/a  ),  Jo  (ka)  ~  1,  and  J\  (ka)  ~  y,  so 


that 


Za  ~ 


2  pc2 
ztua 


Noting  that  pc2  =  /?,  the  adiabatic  bulk  modulus  of  the  fluid  (Kinsler  and  Frey, 
1982;  Junger  and  Feit,  1986),  we  may  write  Za  «  23/iuja.  Now  consider  a  section 
of  the  fluid  cylinder  of  length  l  and  area  2nal.  The  mechanical  impedance  of  the 


section  is  2iralZa ,  thus, 


■> M 


An vl(3 
iu  ’ 


c 

u>  <  - 
a 


Thus  for  very  low  frequencies  the  fluid  cylinder  acts  like  a  spring  of  stiffness 
Airlfi  Newtons  per  meter. 

For  values  of  the  argument  ka  >  27 r,  the  Bessel  function  is  approximated  as 
Jn  (ka)  ~  cos  (ka  —  ^  ~  f )  •  Using  this  relationship,  the  expression  for 

the  spectral  impedance  (with  kz  =  0)  at  high  frequencies  is 


Za  pa  —ipc  cot  ^ ka  —  — ^  ,  ka  >  2n 


The  exact  value  of  normalized  spectral  impedance  from  Eq.  (2.15  )  is  shown 
in  Figure  2.4.  As  in  Case  1,  it  is  purely  reactive,  alternating  between  stiffness  and 
mass  type  reactance.  The  impedance  is  zero  at  the  zeros  of  Jo  (ka) ,  and  infinite 
at  the  zeros  of  J\  (ka). 
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Case  3  -  Ay  Imaginary  When  the  axial  wavenumber  kz  exceeds  the  acoustic 
wavenumber  k,  the  radial  wavenumber  is  imaginary,  and  kr  =  i^Jk2  —  k2  =  ikr. 
For  an  imaginary  argument,  the  Bessel  function  of  the  first  kind  is  replaced  with 
the  modified  Bessel  function  of  the  first  kind.  Thus  Jn  (ikra )  =  inIn  (kra).  Sub¬ 
stituting  in  (2.14),  the  spectral  impedance  is  given  by 


ipck  Iq  (kraj 

Za  —  — - jn  \ 

kr  I\  ( kra) 

For  small  values  of  the  argument,  In  (kTaj 
impedance  is  approximated  as 


(2.16) 


(kra)n  ,  _ 

and  the  spectral 


Za  ~ 


2 18 


ua 


m -i 


Now  when  kz  k,  kT  rs  kz  and  for  large  values  of  kra  the  large  argument 
approximations  for  In  (kraj  yield 


Za 


— pck 
ikz 


ipu 

kz 


(2.17) 


2.2.4  Free  Wave  Propagation 

"We  now  consider  the  propagation  of  free  waves  in  a  fluid  cylinder  for  the  limiting 
cases  of  a  rigid  and  a  free  boundary.  In  contrast  to  the  preceding  discussion 
of  spectral  impedance  in  which  both  the  frequency  and  the  axial  wavenumber 
were  permitted  to  be  free  variables,  we  now  seek  values  of  the  radial  wavenumber 
kr  which  will  satisfy  specific  boundary  conditions.  These  particular  values  of  kr 
determine  the  modes  of  free  wave  propagation.  Again,  the  analysis  is  restricted  to 


the  axisymmetric  case. 
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Free  Radial  Boundary  The  pressure  must  vanish  at  a  free  boundary,  so  we  set 
Pa  =  0  in  Eq.  (2.14).  To  satisfy  this  condition  requres  that  Jo  (kra)  =  0.  Des¬ 
ignating  the  roots  of  this  equation  as  jojn .  the  allowable  values  of  kr  are  given 

by 

^r,m^  jo, mi  “Wl  1, 2,  3,  .  .  . 
and  j0,i  =  2.40,  j0, 2  =  5.52,  j0, 3  =  8.65,  etc. 

For  values  of  m  greater  than  2  the  roots  may  be  approximated  by  j0,m  « 
|(4 m  —  1)  with  an  error  of  less  than  0.1%. 

Each  root  determines  the  characteristics  of  wave  propagation  of  a  particular 
mode  in  the  fluid  cylinder.  Since  the  axial  wavenumber  must  be  real  for  propaga¬ 
tion  to  occur,  each  mode  will  exhibit  a  cutoff  frequency  uJm,  determined  from  the 
dispersion  equation  by  setting  kZym  =  0,  as 

J0,7TIC 
=  — - 

a 

The  phase  velocity  of  the  wave,  c^,  is  computed  from  the  axial  wavenumber 


and  the  group  velocity  cg  is  given  by 

-  —  -fi  1 

“'m  -  dkz,m  ~  \ua  )  J 

The  variation  with  frequency  of  the  phase  and  group  velocities  for  a  fluid 
cylinder  with  a  free  radial  boundary  are  shown  in  Figure  2.5.  This  is  an  example  of 
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’’normal”  dispersion,  in  which  the  phase  velocity  decreases  with  frequency.  (Elmore 
and  Heald,  1985) 

Substitution  of  the  modal  radial  wavenumber  Ay>m  =  jo,m/a  into  the  wave  equa¬ 
tion  yields  the  characteristic  function  describing  the  radial  variation  in  pressure 
for  the  m  th  mode: 

Pm  (p)  =  -Am  Jo  (^jo.m  ~ ^ 


where  Am  is  a  modal  amplitude  coefficient.  The  normalized  modal  amplitude, 
Pm  Mm  is  plotted  as  a  function  of  the  normalized  radius  r/a  in  Figure  2.6  for  the 
first  three  modes.  We  note  that  there  is  no  plane  wave  mode  (m  —  0)  for  this 
case,  since  the  pressure  must  be  zero  at  r  =  a.  For  any  mode,  the  pressure  is  a 
maximum  at  the  center  of  the  cylinder  and  zero  at  the  boundary. 

As  discussed  in  section  1.6,  the  wave  in  the  fluid  cylinder  may  be  constructed 
from  an  infinite  number  of  plane  waves  whose  wavevector  k  makes  an  angle  <f>  with 
the  z  axis.  This  angle  is  determined  from  tan  <f>  =  jr  and  the  angle  for  a  particular 
mode  is 


( fim  =  arctan 


1 

2 


Rigid  Radial  Boundary  The  presence  of  a  rigid  boundary  requires  that  the  radial 
component  of  velocity  be  zero  at  r  =  a.  From  Eq.(2.13)  this  condition  is  satisfied 
when  Ji  ( kra )  =  0.  Letting  jijTn  represent  the  roots  of  this  equation,  modes  exist 
for  values  of  Ay  given  by 


kr,m ®  —  jl,m , 


m  —  0,  1,  2,  . . . 
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where  jlfl  =  0,  =  3.83,  jh2  =  7.01,  ji,3  =  10.17  For  m  greater  than  5 

the  approximation  ji?m  ^  (4m  +  1)  may  be  used  with  an  error  of  less  than  0.1%. 

The  phase  and  group  velocities  for  a  fluid  cylinder  with  a  rigid  boundary  are 
presented  in  Figure  2.7.  The  ji?o  —  0  root  characterizes  the  plane  wave  mode  which 
is  supported  by  the  rigid  wall  boundary  condition.  There  is  no  cut-off  frequency 
for  this  mode,  and  the  phase  and  group  velocities  are  independent  of  frequency 
and  equal  to  the  sound  speed  in  the  fluid. 

The  modeshapes  for  this  case  are  given  by  pm  (r)  =  Am  Jo  (ji and  are 
shown  in  Figure  2.8.  For  the  plane  wave  mode,  the  pressure  is  independent  of  the 
radius. 


2.3  Fluid  Tube 


2.3.1  Introduction 

We  now  consider  a  cylindrical  tube  of  fluid,  infinite  in  length,  with  inner  radius 
a,  outer  radius  fr,  sound  speed  c,  and  density  p.  Such  a  fluid  cylinder  could  be 
physically  realized,  for  example,  by  a  fluid  contained  in  the  annular  space  between 
two  concentric  cylindrical  shells. 

We  seek  to  develop  a  transfer  function  relating  the  state  variables  of  pressure 
P  and  radial  velocity  V  at  the  inner  boundary,  r  =  a,  to  those  at  the  outer 
boundary,  r  =  6,  for  an  arbitrary  excitation  comprising  an  axisymmetric  traveling 
pressure  wave  at  either  boundary.  The  analysis  is  restricted  to  the  axisymmetric 
case  for  clarity,  but  it  may  easily  be  extended  to  the  more  general  non-axisymmetric 
case  using  the  transform  pair  defined  in  equations  2.9  and  2.11  to  describe  the 
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circumferential  variation.  This  will  result  in  an  infinite  sum  of  transfer  functions 
to  represent  the  fluid  tube. 

Figure  2.9  describes  the  model  in  block  diagram  form.  Pa  and  Va  represent  the 
state  variables  of  pressure  and  velocity  at  the  inner  radius,  and  Pb  and  Vb  those  at 
the  outer  radius. 

The  desired  transfer  function  is  a  matrix  of  the  form 


pb 

A 

B 

Pa 

vb 

C 

D 

Va 

where  the  parameters  A,  B,  C.  and  D  are  defined  as: 

A  =  §  k=o>  B  =  ^  | pa=0,  C  =  I va=o  and  D  =  ^  |Pa=0 


2.3.2  Development  of  Transfer  Function 

The  pressure  within  the  fluid  tube  must  satisfy  the  scalar  Helmholtz  equation, 


(Vr  +  k2^  p(kz,<jj,r)  =  0 


where  p(kz  ,u>,r)  is  the  double  Fourier  transform  of  the  pressure  as  described 


in  equation  2.8  and  =  -4-  r-Jr 


dr 


Now  we  assume  a  solution  of  the  form 


p{kz,u,r )  =  [E  40  (krr)  +  F  (krr) 


0i(kzz—ut) 


(2.18) 


where  H\p  ( krr )  is  the  Hankel  function  of  the  first  kind,  order  0,  representing 
waves  traveling  outward  in  the  positive  r  direction,  and  ( kTr )  is  the  Han¬ 
kel  function  of  the  second  kind,  representing  waves  traveling  inward.  A  radial 
wavenumber  kr  is  computed  from  k%  =  k2  —  k2  =  —  fcf ;  and  E  and  F  are 
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complex  coefficients  to  be  determined.  The  radial  velocity  vr  is  found  using  the 
momentum  equation, 


-Vrp(K,u,r)  =  p—vr(kz,u,r)  =  -ipuvr  (kz,u,r) 


Applying  this  to  equation  2.18  results  in 

vr  (kz,cu,  r)  =  —  \E  ( krr )  +  F  H[ 2)  (krr) ] 

uip  L  J 

Defining  the  radial  spectral  impedance  of  the  fluid  in  the  tube  as 


the  radial  velocity  may  be  expressed  as 


vr  (kz,u,r)  =  ±[e  H?  0 krr )  +  F  H (2)  (krv)]  (2.19) 

The  transformed  state  variables  Pa,  Va ,  Pb,  and  14  are  found  by  substituting 
r  =  a  and  r  —  b  into  equations  2.18  and  2.19.  This  results  in  the  following  four 
equations  which  must  be  solved  to  determine  the  complex  coefficients  E  and  F. 
The  exponential  function  exp  i(kzz—UJt)  will  be  suppressed  in  subsequent  equations 
for  brevity. 

Pa  =  E  ( kra )  +  F  H™  ( kra ) 

Va  =  ^[E  (kra)  +  F  H{2)  (kra) 

Pb  =  E  (krb)  +  F  #£2)  (kTb) 

Vb  =  ^\E  (krb)  +  F  H[2)  (Kb) 

Zj 
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These  equations  may  be  written  in  matrix  form  as 


Pb 

1 

E 

— 

Mb 

vb 

L  J 

F 

(2.20) 


and 


Pa 

r  i 

E 

=  Ma 

Va 

F 

(2.21) 


where 


Mn 


and 


Mb 


-| 

t-H 

1 

to 

i 

a\ 

a2 

- 

^H[1]  (kra) 

(kra) 

H(01]  0 Kb ) 

i 

^4° 

h 

K 

* 

FH?  (krb) 

iH?  (krb) 

b3 

K 

(2.22) 


(2.23) 


Multiplying  2.20  by  the  inverse  of  2.23  and  2.21  by  the  inverse  of  2.22  leads  to  the 
elimination  of  the  complex  coefficients  E  and  F,  viz. 

1  -1  f  A  1 

(2.24) 


,  the 


E 

r  -| 

-1 

Pb 

r  i 

-1 

Pa 

= 

Mb 

= 

Ma 

F 

L  J 

.  ^  . 

L 

K 

Multiplying  2.24  by 
identity  matrix,  result  in: 


Mh 


and  noting  that 


r  i 

“ 

-1 

1 

■l" 

Or*1 

1 _ 

Mb 

= 

i 

1 _ 

■ 

A 

= 

Mb 

Ma 

-1 

Pa 

_  A  _ 

va 

By  a  similar  method  we  obtain 


Pa 

1 

r  -] 

~1 

Pb 

= 

Ma 

Mb 

va 

vb 
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Thus  Mb  Ma  is  the  transfer  matrix  from  the  inner  surface  of  the  fluid 
LJL  -]  f  1 — 1 

tube  to  the  outer  surface,  and  M  Mh  is  the  transfer  matrix  from  the 
outer  surface  to  the  inner  surface.  Substituting  the  expressions  for  Ma  and 
Mb  from  (2.22)  and  (2.23)  we  have 

(6jCl4  —  ^2^3)  (~b\(l2  +  &2fll) 

(63&1  —  64O3)  (—63^2  +  b^a\) 

A  a 

where  A  a  —  ai02  —  a2a  and  [Mab]  represents  the  transfer  matrix  from  the  inner 
surface  of  the  fluid  tube  at  r  =  a  to  the  outer  surface,  r  —  b. 

Substituting  the  Hankel  functions  into  the  expression  for  A  a  we  have 

Aa  =  ^  \H i1}  ( kra )  H{2)  ( kra )  -  ( kra )  ( kraj 

The  Wronskian  for  the  Hankel  functions  is 

W  {tf «  (x) ,  H™  (x)}  =  —  =  H%  (x)  H{2)  (x)  -  H™  (x)  (x) 

therefore  Aa  =  and  the  desired  Ay  B ,  C,  D  parameters  are: 

A  =  -ijkra  [h™  ( krb )  H{2)  {kra)  -  ( krb )  H?  (M)]  (2-25) 

B  =  ^ kraZ  [H{Z  (krb)  (kra)  -  ( krb )  (fcra)]  (2.26) 

C  =  7^  \H !(1)  {krb)  h[2)  ( kra )  -  H{2)  {krb)  Ha(1)  (fcra)l  (2.27) 

4  Z  L 

D  =  -t^Ava  [tf{2)  (Av6)  Hf }  (M)  -  #i(1)  (krb)  (*ra)]  (2.28) 


30 


ABCD  Parameters  Expressed  in  Bessel  Functions 
Using  the  definitions  of  the  Hankel  functions,  (x)  =  Jn  ( x)+iYn  (#)  and  H ^  ( x ) 
Jn  (x)  —  %Yn  (. x )  the  ABCD  parameters  may  be  expressed  in  the  form 


A  =  | kra  [Jj  (kra)  y0  (M)  -  J0  (krb)  Yx  (kra)} 

(2.29) 

7T 

B  =  i—kTaZ  [Jo  (fcra)  lo  (&v&)  —  Jo  (fcr6)  Yo  (. kra )] 

(2.30) 

7 T  Jc  d 

c  =  [Ji  (M  n  (Av6)  -  Ji  (M)  H  (M] 

(2.31) 

D  =  T- kra  [  Jx  (Av6)  r0  (M  -  Jo  (M  Yl  (M)} 

Zd 

(2.32) 

Note  that  A  and  D  are  real,  whereas  B  and  C  are  imaginary. 

A  Compact  Notation  for  A  B  C  D  Parameters 
Hermann  and  Mirsky  (1956)  defined  a  new  symbol  to  represent  combinations  of 
Bessel  functions  of  the  form  shown  above.  Their  definition  is: 

Qmn(^r)  Jm  (&r&)  Yn  (/cr5)  Jn  (/tr6)  Ym  ( krCL ) 

Using  this  definition  facilitates  a  more  compact  expression  of  the  ABCD 
parameters,  viz. 

7T 

A  —  ~kra  Qio  (&r) 

7T 

B  —  Q00  (^r) 

c  =  i~Qu(K) 

7 r 

D  —  ~kra  Qoi  (kr) 
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Zero  Thickness  Fluid  Tube 

As  a  check  on  the  preceding  results,  let  a  =  b  ,  that  is,  the  tube  has  zero  thickness. 
Substituting  a  =  b  in  Equation  2.29  we  have 

7T 

A  =  —kTa  [Ji  ( kra )  Y0  ( kTa )  —  Jo  ( kra )  Y\  (A;ra)] 

The  expression  in  brackets  is  the  Wronskian,  which  equals  thus  for  this 
case,  A  —  1.  The  product  terms  in  B  and  C  cancel,  so  that  B  =  C  =  0.  In  similar 
fashion,  D  —  1.  Therefore  the  transfer  matrix  reduces  to 


A 

1  0 

Pa 

vb 

0  1 

Va 

That  is,  Pb  —  Pa  and  Vj,  =  Va  ,  as  expected. 

A  B  C  D  For  Imaginary  Radial  Wavenumber 
When  the  axial  wavenumber  kz  exceeds  the  acoustic  wavenumber  k.  the  radial 
wavenumber  is  imaginary.  Using  the  same  notation  as  in  Case  3  for  the  fluid 
cylinder  we  substitute  kr  =  ikr  into  Equations  2.25  to  2.28  and  replace  the  Hankel 
functions  with  the  modified  Bessel  functions,  viz. 

(iKa)  =  ^n+1  Kn  ( kTa ) 

and 

2 

H(rp  ( ikra )  =  2 inIn  (; kra )  -  —p^Kn  (kra) 

The  resulting  ABC  D  values  for  the  fluid  tube  with  imaginary  radial  wavenum¬ 
ber  are 

Ko  (krb)  h  (kra)  +  K0  (kra)  h  (m)] 


A  —  kra 


(2.33) 
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B  =  ipckra  [Ko  (kraj  I0  (krbj  —  Kq  (^krb)  Iq  (kTaj  (2.34) 

■t2  _  _  - 

C  =  t-^[K1  (krb)  h  ( %a )  -  Kx  (kra)  h  (m)]  (2.35) 

D  =  kTa  K0  (kraj  I\  (krbj  +  K\  (krbj  I0  (kraj^  (2.36) 

Again  we  note  that  A  and  D  are  real,  whereas  B  and  C  are  imaginary. 

2.3.3  Transfer  Function  From  b  to  a 

As  mentioned  earlier,  the  transfer  function  from  the  outer  surface  at  r  =  b  to  the 
inner  surface  at  r  =  a  is  defined  by 

[Mba\  =  [Ma]  [M,]-1 

When  this  is  solved  for  the  A  B  C  D  parameters  in  the  same  manner  as  for 
the  preceding  case,  the  expressions  take  the  same  form  but  with  the  values  of  a 
and  b  interchanged,  viz. 

A  =  | krb  [Jj  (krb)  V0  (kra)  -  J0  (kra)  Y1  ( krb )]  (2.37) 

B  =  i—krbZ  [Jo  ( kTb )  Yo  ( kra )  —  Jo  (kTa)  Yo  (Kb)]  (2.38) 

C  =  [Jx  (krb)  Yx  (kra)  -  Jx  (kra)  Yx  ( krb )]  (2.39) 

Z  Zj 

D  =  ^krb  [Jr  (kra)  Y0  (krb)  -  Jo  (krb)  Yi  (kra)] 

& 


(2.40) 
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2.3.4  Large  Argument  Approximations  for  Transfer  Functions 

The  A  B  C  D  parameters  for  large  values  of  the  argument  are  determined  using 

the  asymptotic  values  of  the  Hankel  functions, 


(z) 


•  /  _  ^7 r  7T  \ 

,i[z- nj-lj 


Z  »  27T 


and 


(z) 


i(z-nf-f) 


Substituting  these  relationships  into  Equations  2.25  to2.28  and  simplifying  the 
results  yields  the  following  expressions  for  A  B  C  D  for  large  real  values  of  the 
non-dimensional  radial  wavenumbers,  kra  and  krb  : 


A 


cos  kr  ( b  —  o) 


(2.41) 


B 


(2.42) 


C 


i  fa 

ZVb 


sin  kr  ( b  —  a) 


(2.43) 


D 


cos  kr  ( b  —  a) 


(2.44) 


Now  if  we  let  h  represent  the  thickness  of  the  fluid  tube,  that  is  h  =  b  —  a,  and 
assume  that  h  is  very  small  compared  to  b,  so  that 
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then 

A  —  D  rs  cos  krh 

B  ~  iZ  sin  krh 

C  R*  —  sin  krh 
Zt 

Recalling  that  Z  —  k  =  and  kr  =  \Jk2  —  k%,  we  observe  that  the 

asymptotic  expressions  for  the  two  port  model  of  a  fluid  tube  reduce  to  those 
developed  for  plane  waves  in  another  study.  (Stepanishen,  1985) 


Imaginary  Large  Argument 

When  the  non-dimensional  radial  wavenumber  is  large  and  imaginary,  kra~A>  2ir. 

exp  (for  ft) 


the  modified  Bessel  functions  are  replaced  with  In(kra ) 

vE  exp 


.  and  Kn  (kra) 

\/  2irfcrO  ' 


Using  these  expressions  in  Equations  2.33  to  2.36  yields  the  following  equations 
for  the  A  B  C  D  parameters  for  large  imaginary  values  of  the  argument  kra  : 

A  =  D  r*  a/~~  cosh  —  l^j  exp  ^/craj 


_  ipcA:  [a  b 
B  r*  -=r—  i/t  smh( - l)exp(fcra) 


kr 


a 


sinhf--l)exp 
pck  Vo  \  a  / 


We  note  that  these  expressions  indicate  that  the  negative  square  root  must  be 
chosen  for  kr  in  order  to  represent  an  evanescent  wave  which  decays  with  radial 


distance. 
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2.3.5  Discussion  of  Fluid  Tube  Transfer  Function 

We  now  consider  the  behavior  of  the  fluid  tube  A  B  C  D  parameters  as  functions 
of  the  radial  wavenumber  and  the  fluid  tube  radial  dimensions. 

It  is  convenient  to  begin  with  the  asymptotic  approximations  for  large  ar¬ 
guments  since  the  trigonometric  functions  involved  are  simpler  than  the  exact 
relationships  which  are  the  products  of  Bessel  functions. 

We  also  introduce  the  dimensionless  ratio  of  the  outer  to  inner  radius,  s  =  b/ a. 
Substituting  into  Equations  2.41  to  2.44  we  have 

1 

A  =  D  ~  —=  cos  [kra  (s  —  1)] 

Vs 


B 


ipck 

kry/s 


sin  [kra  (s  —  1)] 


C 


ikr 

pcky/s 


sin 


[kra  (s  —  1)] 


These  functions  are  periodic  with  period  krtns  (s  —  1)  =  27m.  The  amplitude 
coefficient  1/y/s  =  y/§  reflects  the  decrease  in  magnitude  of  the  pressure  and 
velocity  at  b  relative  to  a  due  to  cylindrical  spreading. 

As  s  is  increased  the  frequency  of  variation  with  radial  wavenumber  kra  also 
increases;  that  is,  there  is  a  larger  number  of  pressure  and  velocity  cycles  between 
the  inner  and  outer  radii. 

Now  let  the  axial  wavenumber  kz  =  0  so  that  there  is  no  propagation  in 
the  axial  direction.  For  this  case  the  radial  wavenumber  kT  equals  the  acoustic 
wavenumber  k  .  The  equations  of  state  for  this  case  are: 
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Pb  =  —%  cos  [ka  ( s  —  1)]  +  i^y=Va  sin  [ka  (s  —  1)] 
v  s  \/s 


Vb  —  —  sin  [ka  (s  —  1)1  H — p  cos  [ka  (s  —  1)1 
pc^/s  y/S 

Now  when  ka  (s  —  1)  =  m r,  n  =  0, 1, 2, ...  the  sine  terms  vanish,  and  the 

cosine  terms  equal  (— l)n  .  The  transfer  equations  reduce  to 


A  =  (-1)' 


Vt  =  (-l)”l£= 

Vs 

Thus  the  pressure  and  velocity  at  b  are  either  in  phase,  or  180°  out  of  phase 
with  those  at  a.  This  condition  occurs  when  the  tube  thickness,  h  =  (6  —  a)  is  an 
integral  number  of  half  wavelengths,  h  =  n\n/2  or  at  a  frequency  un  =  Tmcjh. 
Similarly,  when  ka  (s  —  1)  =  (2 n  +  1)  7t/2,  the  sine  terms  equal  (— l)n  and  the 


cosine  terms  vanish.  The  result  is 


Pb  =  i(-l)npc- 


vb  =  i(-iy 


pCy/s 


This  will  be  true  when  the  thickness  h  is  an  odd  number  of  quarter  wavelengths, 
i.  e.  h  =  (2 n  +  1)  A/4  ,  or  un  =  (2 n  +  1)  7tc/2/i 


The  conditions  just  discussed  describe  a  radial  standing  wave  in  the  fluid  tube. 
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For  the  general  case  of  interest  in  which  the  large  argument  approximations 
are  not  valid,  the  A  B  C  D  parameters  are  computed  using  one  of  the  sets  of 
expressions  presented  earlier,  such  as  Equations  2.29  to  2.32.  The  parameters  are 
plotted  in  Figure  2-10  as  a  function  of  the  non-dimensional  wavenumber  kra  for 
three  values  of  the  ratio  of  radii,  s  =  1.5,  2,  and  3.  The  plots  for  B  and  C  are 
normalized  with  respect  to  the  spectral  impedance  Z,  that  is,  C'  —  —iCZ  and 
B'  =  - iB/Z . 

The  behavior  of  the  A  B  C  D  parameters  in  gross  aspect  follows  the  pattern 
indicated  by  the  asymptotic  approximations.  Both  A  and  D  are  similar  to  cosine 
functions,  whereas  B  and  C  resemble  sine  functions.  The  amplitude  of  the  envelope 
of  the  parameters  decreases  with  increasing  values  of  s. 

Using  the  asymptotic  expressions  for  the  Bessel  functions  with  small  argument 
in  Equations  2.29  to  2.32  and  taking  the  limit  as  kr  — *  0,  we  find: 

lim  A  —  1:  lim  B  =  lim  C  =  0;  lim  D  —  - 

kr-^O  kr-*  0  kr-*  0  fcr->0  S 

This  condition  will  occur  when  the  axial  wavenumber  kz  equals  the  acoustic 
wavenumber  k.  For  this  case,  the  pressures  at  a  and  b  are  equal,  and  V&  =  VaJs. 
This  would  be  the  situation  for  a  plane  wave  propagating  in  the  fluid  tube  with 
no  radial  component. 

Figures  2-11  to  2-14  depict  the  A  B  C  D  parameters  as  three  dimensional 
surfaces  for  values  of  kra  ranging  from  0  to  5,  and  for  s  ranging  from  2  to  5.  The 
features  described  above  are  clearly  indicated  in  these  plots. 
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2.3.6  Free  Wave  Propagation 

We  now  consider  the  propagation  of  free  waves  in  a  fluid  tube.  In  contrast  to  the 
preceding  analysis,  in  which  both  the  the  axial  wavenumber  kz  and  frequency  were 
assumed  to  be  independent  variables,  now  their  relationship  is  determined  by  the 
radial  wavenumber  kr,  which  in  turn  takes  on  discrete  values  that  are  established 
by  the  boundary  conditions. 

The  transfer  function  parameters  provide  a  convenient  means  for  determining 
the  conditions  which  govern  the  propagation  of  free  waves  in  the  fluid  tube.  We 
recall  that  the  transfer  equations  are: 


Pb  —  PaA  +  VaB  (2.45) 

Vb  =  PaC  +  VaD  (2.46) 

in  which  Pa.  Pb)  Va  and  Vb  are  the  state  variables  at  the  surfaces  r  =  a  and 
r  =  b,  respectively.  Four  cases  for  free  wave  propagation  may  be  defined  in  terms 
of  the  state  variables  at  the  boundaries: 


Case 

Description 

Boundary  Conditions 

I 

Free-Free 

II 

o 

II 

o 

II 

Free-Fixed 

pa  =  0,  Vb  =  0 

III 

Fixed-Free 

Va  =  0,  Pb  =  0 

IV 

Fixed-Fixed 

II 

o 

& 

II 

O 
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Case  I,  Free-Free  Boundary 

Now  consider  Case  I.  For  the  free-free  case,  the  pressure  must  be  zero  at  both  the 
inner  and  outer  surfaces  of  the  fluid  cylinder.  Setting  Pa  =  Pb  =  0  in  Equations 
2.45  and  2.46  results  in  the  following  two  equations:  VaB  =  0  and  Vb  =  VaD. 
This  in  turn  requires  that  B  =  0.  Substituting  for  B  in  Equation  2.26  yields  the 
characteristic  equation  for  the  propagation  of  free  waves: 


( kras )  i42)  M  -  #0°  ( has )  (kTa)  =  0  (2.47) 

Denoting  the  roots  of  this  equation  as  hm,  free  waves  may  exist  for  values  of 
the  radial  wavenumber  such  that 

,  _  hm 

K>rm 

a 

To  investigate  these  roots,  we  first  consider  the  case  for  which  the  argument 
kr,ma  is  large.  Using  the  asymptotic  expressions  for  the  Hankel  functions,  the  roots 
of  Equation  2.47  occur  at  values  of  kr  which  satisfy  the  following  relation: 


exp  [i  (2 kra  {5  —  1})]  =  1 


This  condition  is  met  for  values  of  the  radial  wavenumber  kr  m  such  that 


7  rrm 

k>r,m  —  7  77 

a  (s  —  1) 


,  m  =  1,2,3,  ..  . 


Expressing  this  in  terms  of  the  tube  thickness  (b  —  a) , 


k  — 

*Xjr,m 


rrm 

{b-a) 


(2.48) 
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we  see  that  the  radial  wavenumber  kr^m  for  a  given  mode  number  m  decreases 
as  the  tube  thickness  increases. 

The  exact  values  of  the  roots  hm  must  be  determined  by  numerical  solution 
of  the  characteristic  equation,  2.47.  The  asymptotic  value  of  the  root  kTiTna  for  a 
particular  mode  number  m  is  used  as  the  estimated  root  in  a  root  solving  algorithm. 
The  roots  for  the  first  ten  modes  are  listed  in  Table  I  for  four  values  of  s  :  1.5, 2, 3, 
and  5.  Also  shown  in  Table  I  is  the  error  (expressed  as  a  percentage  of  the  exact 
value)  which  results  when  the  asymptotic  approximation  for  the  root  is  used.  For 
a  particular  mode  number,  the  error  increases  as  the  ratio  s  —  bja  increases. 

The  error  is  plotted  as  a  function  of  the  mode  number  for  four  values  of  s  in 
Figure  2-15. 

Case  II,  Free-Fixed  Boundary 

In  Case  II,  the  inner  surface  of  the  fluid  tube  is  free,  and  the  outer  surface  is 
blocked.  Thus  in  Equations  2.45  and  2.46,  Pa  and  Vb  are  zero.  This  is  satisfied  if 
the  parameter  D  equals  zero.  From  Equation  2.28,  the  characteristic  equation  for 
Case  II  is 

H[2)  ( kras )  h£]  (kra)  -  H a(1)  (kras)  ( kra )  =  0  (2.49) 

Again,  using  the  asymptotic  expressions  for  the  Hankel  functions  to  determine 
the  approximate  roots  of  2.49  we  obtain 

exp[i2(kra{s  -  1}  -  ^)]  =  1 
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Thus  the  approximate  values  of  the  roots  for  this  case  are  given  by 

hm  =  kr  ma  -  (2"  +  11);,  m  =  1,  2,  3,...  (2.50) 

2(s-l) 

This  result  was  used  in  the  characteristic  equation  to  determine  the  exact 
values  of  the  roots  for  the  first  ten  modes  as  in  Case  I.  The  results  are  presented 
in  Table  II.  Figure  2-16  shows  the  error  as  a  function  of  the  mode  number  for  four 
values  of  s. 

Case  III,  Fixed-Free  Boundary  Case  III  is  the  reverse  of  Case  II  in  that  the  inner 
boundary  is  rigid  and  the  outer  boundary  is  free,  so  that  Va  and  are  zero.  Thus 
A  =  0  and  the  characteristic  equation  is 

(Kas)  H[2)  ( kra )  -  H(2)  ( kras )  (kra)  =  0  (2.51) 

For  large  values  of  the  argument  kra  the  asymptotic  form  of  this  equation  is 

exp[i2(kra{s  -  1}  -  -)]  =  1 

This  is  identical  with  the  result  for  Case  II,  hence  the  approximate  values  for 
the  roots  of  Equation  2.51  are  given  by  Equation  2.50.  Table  III  lists  the  exact 
values  of  the  roots  for  the  first  ten  modes  of  Case  III  for  four  values  of  s.  Figure 
2-17  shows  the  error  as  a  function  of  mode  number  for  the  same  values  of  s. 

Comparing  Figure  2-16  with  Figure  2-17  we  note  that  the  approximate  formula 
predicts  a  value  of  the  roots  which  is  too  large  in  Case  II,  and  too  small  in  Case 


III. 
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Case  IV,  Fixed-Fixed  Boundary  Finally,  for  Case  IV,  the  fixed  boundaries  at  a 
and  b  imply  that  the  radial  velocities  Va  and  Vb  vanish,  which  in  turn  requires  that 
C  =  0.  The  characteristic  equation  for  this  case  is 

[. H j(1)  (kras)  H{2)  (kra)  -  H[2)  (kras)  H[1}  (fcra)]  =  0  (2.52) 

pcKd 

The  asymptotic  solution  to  this  is  identical  with  that  for  Case  I,  and  the 
approximate  roots  are  given  by  Equation  2.48.  The  exact  roots  for  the  first  ten 
modes  are  listed  in  Table  IV,  and  Figure  2-18  shows  the  error  as  a  function  of 
mode  number  and  s.  Comparing  the  results  for  Case  I  with  those  for  Case  IV,  we 
note  that  the  asymptotic  formula  yields  a  positive  error  for  the  free-free  case,  and 
a  negative  error  for  the  fixed-fixed  case.  We  note  that  Case  IV  is  the  only  one  to 
which  supports  a  plane-wave  mode. 

Up  to  this  point  in  the  analysis  we  have  used  the  characterictic  equations  in 
complex  form  for  ease  of  manipulation  in  determining  the  asymptotic  expressions 
for  the  roots.  The  AB  C  D  parameters  may  also  be  cast  in  the  form  given  in  Equa¬ 
tions  2.29  to  2.32,  which  results  in  the  following  expressions  for  the  characteristic 
equations: 

Case  I: 

J0  ( kra )  Y0  ( kras )  —  J0  ( kras )  Yq  ( kTa )  =  0 

Case  II: 

Ji  ( kras )  Yq  ( kTa )  —  Jq  ( kra )  Y\  ( kras )  =  0 

Case  III: 

Ji  ( kra )  Yq  ( kras )  —  Jo  ( kras )  Y\  ( kra )  =  0 
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Case  IV: 

^£[Ji  (kTa)  Vi  (kras)  -  J:  (kras)  Vj  (fcro)]  =  0 
pcfca 

The  characteristic  equations  for  the  four  cases  analyzed  for  the  fluid  tube  are 
plotted  as  a  function  of  the  non-dimensional  wavenumber  kra  in  Figures  2-19  to 
2-22  for  radii  ratios  of  1.5,  2,  and  3.  We  note  that  as  the  ratio  of  radii  5  increases, 
the  number  of  zeros  for  a  given  range  of  the  non-dimensional  wavenumber  also 
increases.  Physically  this  simply  means  that  the  modal  density  is  greater  in  a 
thick  walled  fluid  tube  than  in  a  thin  one. 

2.3.7  Mode  Shapes  in  a  Fluid  Tube 

We  now  consider  the  mode  shapes  of  the  pressure  and  velocity  in  a  fluid  tube  with 
the  four  sets  of  boundary  conditions  described  above.  Guidance  for  this  effort  is 
provided  by  the  analysis  of  a  rectangular  fluid  waveguide  (Redwood,  1960).  The 
rectangular  waveguide  is  the  limiting  case  of  a  fluid  tube  of  infinite  radius.  The 
mode  shapes  in  the  rectangular  waveguide  comprise  simple  sinusoidal  functions, 
considerably  simpler  than  those  of  the  finite  radius  fluid  tube. 

To  proceed,  we  use  Equation  2.18  for  pressure  and  Equation  2.19  for  the  radial 
velocity.  Using  the  appropriate  boundary  conditions  in  these  equations  extablishes 
the  relation  between  the  unknown  coefficients  E  and  F 

For  example,  in  Case  I,  substitution  of  r  =  a  and  p  =  0  yields 

F  „  (M 
H,'2'  (h,„) 


where  hm  =  kr^ma  is  the  m  th  root  of  the  characteristic  equation. 
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Thus  the  pressure  for  the  m  th  mode  is  given  by 


pm(r,uj,z)  =  E 


H ^  (hm) 
H f  (hm) 


-Ut)\ 


where  the  axial  and  time  dependence  have  been  included. 

The  Hankel  functions  may  be  written  in  polar  form  as  (hm)  =  Mn  (hm)  exp  (i@m) 
and  Hp  (hm)  =  Mn  (hm)  exp  (— *0m) ,  where  Mn  (hm)  =  ±\/jn  (hm  f  +  Yn  (hm)2 
and  0m  =  arctan  ( '  Then  =  exp  (2 i6m)  ,  and  the  pressure  may  be 

written  as: 


Pm  ( r,u,z )  =  E 


hP  (k.1) 


(2) 


hm-  )  exp  (2 iOm) 


exp  [i  ( kzz  —  cut)]  (2.53) 


Equation  2.53  expresses  the  solution  for  pressure  in  complex  notation.  To  find 
the  actual  physical  pressure,  we  replace  the  Hankel  functions  with  the  Bessel  and 
Neuman  functions,  and  the  complex  exponentials  with  sine  and  cosine  functions. 
Expanding  the  result  and  taking  the  real  part  results  in: 


Pm  (r,U,z)  —  r  |  Jo  (km*-j  [1  -  cos(20m)]  -  Yft  (^™“)  Sin  (28m)  ]  (cos  (kzz  -  U>t)} 

(2.54) 

The  variable  E  is  seen  to  be  a  real,  arbitrary  pressure  amplitude  coefficient. 
The  radial  velocity  is  obtained  using  Equation  2.19, 

vr  (r,  co,  z)  =  E  ( krr )  +  F  ( kTr )  exp  [i  ( kzz  —  cut)] 

pck  L 

Using  the  expression  for  F  found  above,  the  radial  velocity  for  the  m  th  mode 
is  given  by 


ikrE  f 

V""lr-U’’z)  =  -M\ 


H\l)  (^m~)  “  ]  exp  (2Wm)  exp  [?  (kzz  -  \ 
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Expanding  and  taking  the  real  part,  the  final  expression  for  the  radial  velocity 
is 


vr,m  (r,u,z)  =  ^  | h  (hm^j  sin  (20m)  -  Fi  (hm^j  [1  +  cos  (20m)]}  {cos  %z  ~  w*)} 

The  axial  velocity  may  be  obtained  by  applying  the  momentum  equation  to 
Equation  2.54  with  the  result 

k 

vz,m  (r,  U,z)  =  Pm  (r,u,z) 

The  axial  velocity  is  real  and  in  phase  with  the  pressure,  which  one  would 
expect  for  a  propagating  wave. 

The  equations  for  the  mode  shapes  of  the  remaining  three  cases  are  identi¬ 
cal  with  those  for  Case  I,  but  with  the  appropriate  root  from  the  characteristic 
equation,  and  with  the  angle  6m  determined  by  the  boundary  conditions. 

Thus  for  Case  II  (free-fixed),  setting  p  —  0  at  r  =  a  results  in  the  same  value 
of  0rn  as  was  found  for  Case  I. 

The  boundary  condition  for  Case  III  (fixed-free)  requires  that  p  =  0  at  r  =  b 


with  the  result 


6m  =  arctan 


Yq  (hmS) 


L  ^0  (fomS)  J 

In  Case  IV  (fixed-fixed) ,  setting  vr  =  0  at  r  =  a  yields 

'n  (hmy 


dm  =  arctan 


J\  (hm)  _ 

Plots  of  the  pressure  and  radial  velocity  as  a  function  of  the  normalized  radius 
(u  =  for  the  first  three  modes  for  each  of  the  cases  studied  are  presented  in 
Figures  2-23  to  2-26.  The  radius  ratio  s  =  10  for  these  plots. 
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The  modes  in  the  fluid  tube  are  considerably  more  complicated  than  those  in 
a  rectangular  fluid  waveguide  or  in  a  fluid  cylinder.  This  is  primarily  because  the 
boundary  conditions  for  the  fluid  cylinder  are  specified  at  two  different  physical  di¬ 
mensions,  the  inner  and  outer  radii  (or,  equivalently,  the  inner  radius  and  the  radii 
ratio),  whereas  the  rectangular  waveguide  and  the  fluid  cylinder  are  characterized 
by  one  physical  dimension.  In  the  case  of  the  rectangular  waveguide,  it  is  the 
thickness  h,  and  for  the  fluid  cylinder,  the  radius  a.  In  the  rectangular  waveguide, 
the  modes  are  described  by  trigonometric  functions  (sines  and  cosines),  and  the 
peak  modal  amplitudes  are  ±1  for  all  cases.  Also,  there  is  ’’left-right”  symmetry 
between  the  free-fixed  modes  of  Case  II  and  the  fixed-free  modes  of  Case  III.  That 
is,  interchanging  the  boundaries  results  in  a  reflection  of  the  mode  shapes  about 
the  center  of  the  waveguide.  This  is  not  the  case  for  the  fluid  tube. 

Mode  shapes  for  the  free-free  situation  of  Case  I  are  shown  in  Figure  2-23. 
There  is  no  ’’plane  wave”  (m  =  0)  mode,  unless  one  considers  the  trivial  case  for 
which  the  pressure  is  everywhere  zero. 

A  salient  feature  of  these  plots  is  that  the  peak  modal  amplitudes  of  both 
pressure  and  radial  velocity  vary  as  functions  both  of  the  mode  number  and  the 
radii  ratio,  s. 

For  small  values  of  s,  eg  s  =  2,  the  mode  shapes  appear  almost  sinusoidal. 
As  s  increases,  the  mode  shapes  appear  more  like  the  Bessel  functions,  with  large 
peak  amplitudes  near  the  inner  radius,  and  decreasing  in  amplitude  as  the  outer 
radius  is  approached.  With  increasing  s,  the  node  for  velocity  in  mode  1  moves 
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towards  the  inner  radius  a. 

The  fluid  cylinder  with  both  boundaries  fixed  (Case  IV)  supports  a  plane 
wave  mode,  unlike  the  other  three  cases.  This  is  evident  by  inspection  of  the 
characteristic  equation  for  Case  IV,  which  exhibits  a  root  at  kra  =  0.  Mode  shapes 
for  m  =  1  to  m  =  3  are  shown  in  Figure  2-24.  As  s  increases,  the  node  for  the 
pressure  in  mode  1  moves  towards  the  outer  radius  b. 

Cutoff  Frequencies 

All  of  the  modes  for  the  fluid  cylinder  (with  the  exception  of  the  plane  wave  mode 
for  case  IV)  will  exhibit  a  cutoff  frequency  below  which  propagation  will  not  occur. 
This  cutoff  frequency  is  established  by  the  root  of  the  characteristic  equation  for 
the  mode,  hm,  the  inner  radius  a,  and  the  sound  speed  c: 

hmc 

U)  c  — 

a 

2.4  External  Fluid 


2.4.1  Introduction 

The  final  component  in  any  arrangement  of  concentric  cylindrical  models  is  the 
external  medium.  In  contrast  to  the  other  components,  which  have  finite  radial 
boundaries,  the  external  medium  is  presumed  to  extend  to  infinity  in  the  radial 
direction.  This  assumption  permits  the  radiation  of  acoustical  energy  in  the  radial 
as  well  as  in  the  axial  direction. 

In  this  study  we  restrict  our  attention  to  inviscid  fluid  media,  such  as  air  and 
water.  However,  solid  external  media  could  be  incorporated  if  desired  by  using 
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elasticity  theory. 

A  number  of  authors  have  discussed  the  propagation  of  cylindrical  waves,  in¬ 
cluding  Morse  (1948),  Rschevkin  (1963),  Morse  and  Ingard  (1968),  and  Junger  and 
Feit  (1986). 

An  especially  cogent  analysis  by  Junger  (1953)  provides  an  excellent  paradigm 
for  the  present  discussion.  In  his  paper,  Junger  postulates  an  infinitely  long  cylin¬ 
drical  surface  which  experiences  a  time  harmonic,  spatially  periodic  standing  wave 
excitation  in  both  the  circumferential  and  axial  directions.  He  then  constructs 
solutions  to  the  inhomogeneous  wave  equation  and  obtains  expressions  for  the 
specific  acoustic  impedance  in  the  external  fluid.  In  a  section  titled  ’’Physical  In¬ 
terpretation”  Junger  presents  an  especially  lucid  description  of  cylindrical  waves 
in  the  fluid  medium. 

In  this  section,  we  analyze  the  propagation  of  cylindrical  waves  in  the  external 
medium.  The  pressure  and  particle  velocity  are  described  at  the  medium  interface 
and  in  the  field  as  functions  of  frequency  (acoustic  wavenumber)  and  the  axial 
wavenumber,  kz.  The  spectral  impedance  at  the  interface  is  also  discussed. 

2.4.2  Physical  Description 

Consider  an  infinite  length  cylinder  of  radius  a  immersed  in  an  inviscid,  homo¬ 
geneous  fluid  which  extends  radially  to  infinity.  A  travelling  harmonic  pressure 
wave  of  amplitude  Pq,  frequency  a?,  and  axial  wavenumber  kz  is  postulated  at  the 
surface  of  the  cylinder,  as  follows: 


pn(a,  0,z.t)  =  P0  cos{n6)ei('kzZ~u}t) 
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The  waves  in  the  external  fluid  medium  must  decay  to  zero  at  r  =  oo  and  also 
must  satisfy  the  reduced  Helmholtz  equation, 


(V2  +  k2)p(r,e,z,t)=  0 


w 


here 


9  d 2  Id  Id 2  d 2 

V2  — - 1 - 1 - 1 - 

dr 2  r  dr  r2  d62  dz 2 


Assume  a  field  solution  of  the  form 


Pn(r,  6,  z,  t )  =  cos  n6  AnH^  ( kTr )  +  BnH £2)  (kTr)  ev 


0i(kzz—ujt) 


For  large  values  of  the  argument  kTr  this  may  be  written  as 


pn(r.0.zA)  =  cos(ne)J—-\Anei{krr+kzZ  ™  +  Bnei(kzZ  Kr 

V  7TKrT  *• 

The  first  term  within  the  brackets  describes  a  wave  propagating  in  the  positive 
axial  and  radial  directions,  that  is,  an  outgoing  conical  wave.  The  second  term 
represents  an  incoming  wave.  Since  we  are  interested  only  in  the  outgoing  wave, 
we  discard  the  second  term,  and  the  Bn  coefficients  are  set  equal  to  zero.  Thus 
the  pressure  in  the  field  is  given  by 


pn(r,  0,  z,  t )  =  cos  (riff)  An  H^r>  ( krr )  et(-kzZ 


(2.55) 


Equating  this  expression  to  the  excitation  pressure  wave  at  r  =  a  enables 
determination  of  the  coefficients  An  as: 


An  — 


Po 


hL1]  (i kra ) 
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The  velocity  in  the  field  is  determined  using  the  harmonic  form  of  the  momen¬ 
tum  equation, 


Vp  =  ipckv 


In  cylindrical  coordinates  this  is  written 

v=_Lfe^  +  5^  +  e&" 

ipck  [  Or  r  m  *dz_ 

The  three  components  of  velocity  thus  are: 
Radial: 


=  ervr  +  eeve  +  ezvz 


vr  =  - 


ikr  cos  (n#)  Po 
Ki) 


pckHT  ( kra )  L  kTr 


n 


Hn}  (krr)  ~  H(n+ 1)  (M 


0i(kzz—u)t) 


(2.56) 


Tangential: 


Axial: 


„  insm(n9)P0HW  (Ayr)  i{kzZ^t) 

Uq  /-x  c 

pckrHk  '  ( kra ) 


h  cos (n0)PoH^  (M_eiikzZ^t) 
pckH^  (kra) 


2.4.3  Wave  Propagation 


Real  Radial  Wavenumber  For  clarity  we  consider  the  axi symmetric  case,  n  =  0, 
to  investigate  the  propagation  of  cylindrical  waves  in  the  external  medium. 

Thus  for  n  =  0  the  pressure  in  the  field  is  given  by: 

p(r,z,t)  - 

#o  (M) 


51 


Since  n  =  0  the  tangential  particle  velocity  is  zero.  The  radial  particle  velocity 


is  given  by: 


_  ikrP0H[l)  (Ayr)  AkzZ_ 


Vr  = 


c ot) 


pckH, ( Aya ) 
and  the  axial  particle  velocity  is 

_  ^zPqHq  ^  (Ayr)  Mkzz-ut) 

pckH^  ( kra ) 

The  Hankel  functions  may  be  written  in  polar  form  as 

Hi1]  (z)  =  Mn(z)ei0n(z) 

where  Mn(z)  =  \J Jn  (zf  +  Yn  (zf  and  9n{z)  =  arctan  .  With  this  nota¬ 


tion, 


p(r,z,t)  = 


PoMo(A:rr)ei(0o(fcrr)“0o(Ava)+fe2_“t) 


M0(kra ) 

For  large  values  of  the  argument  kTr ,  90  (krr)  — >  Ayr  —  |  and 


Mo(krr )  — »  thus  the  pressure  in  the  hmit  becomes 

p(r,  z ,  t)  - 


Po 


V  7tkrr  M0(kra ) 

where  <j>  —  60  (kra  +  is  a  fixed  phase  angle  This  describes  a  traveling  wave 
with  wavevector  components  Ay  and  Ay .  The  direction  of  propagation  is  indicated  by 
the  wavevector  k  =  er  Ay+e,Ay  which  is  inclined  at  an  angle  Q-k  —  arctan  v~  with 
respect  to  the  z  axis. 

The  wavefront  (surfaces  of  constant  phase)  are  conical,  symmetrical  with  re¬ 
spect  to  the  z  axis,  and  perpendicular  to  the  wavevector  k.  These  relationships 


are  shown  in  Figure  2-25. 


If  the  axial  wavenumber  kz  is  zero,  the  radial  wavenumber  kr  equals  the  acoustic 
wavenumber  k.  This  represents  a  cylindrical  wave  traveling  in  the  radial  direction, 


lim  p{r,t) 

KrT—^CO 


Po 


nkrr  M0(kra ) 


^(hrT— U)t—<j>) 


Note  that  the  amplitude  decays  as  which  is  characteristic  of  cylindrical 
waves.  This  also  satisfies  the  requirement  that  the  wave  vanish  at  r  =  co.  We 
may  compare  this,  for  instance,  with  a  spherical  wave  which  falls  off  as  K  The 
wavefronts  are  cylinders  concentric  with  the  z  axis. 


Zero  Radial  Wavenumber  As  the  axial  wavenumber  approaches  the  acoustic 
wavenumber  k  the  radial  wavenumber  approaches  zero  and  radial  phase  veloc- 
ity,  jr,  becomes  infinite.  The  result  is  a  purely  plane  wave  propagating  in  the  z 
direction. 


Imaginary  Radial  Wavenumber  When  the  axial  wavenumber  exceeds  the  acoustic 
wavenumber,  the  radial  wavenumber  takes  on  a  purely  imaginary  value.  We  will  use 
the  notation  kT  to  represent  the  magnitude  of  the  imaginary  radial  wavenumber, 
that  is, 

kr  =  ikr  =  i\jk\  —  k2,  kz>  k 

Using  the  expression  for  Hankel  functions  of  imaginary  argument,  H W  (ix)  = 
Z^-(n+i) Kn  _  the  expression  for  the  pressure  is 
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For  large  arguments,  Kn  (x)  — ►  ^e~x  thus  for  krr  large,  the  pressure  becomes 


Po 


lim  p(r,z,t )  =  x  ii  nl 

krr-*ocyK  K0(kra)\2hr 


which  represents  an  evanescent  wave  propagating  in  the  z  direction  and  falling 
off  rapidly  in  the  radial  direction  as  e~krT . 


2.4.4  Radial  Spectral  Impedance  of  External  Fluid 

We  wish  to  study  the  radial  spectral  impedance  of  the  external  fluid  as  a  function  of 
the  axial  wavenumber  kr  and  the  acoustic  wavenumber  fc.  The  general  expression 
for  the  radial  spectral  impedance  of  the  external  fluid  may  be  obtained  by  dividing 
equation  2.55  by  equation  2.56,  with  the  result: 


Zryn 


ipckH W  (kTa) 


K  [£-rHP  (kr r)  -  H$+l)  (krr) 

Specializing  this  to  the  n  =  0  case  and  evaluating  at  r  =  a  results  in 


Za  =  ~~ipC 


.i-I 


(  kz^ 

f, 


4°  (Ka) 

H?  (Ka) 


(2.57) 


The  nature  of  Za  thus  depends  on  the  values  of  kz  and  k.  Three  cases  of 
interest  include: 

(I)  kz  =  0;  kr  —  k  (cylindrical  wave) 

(II)  kz  <  k;  kr  is  real  (high  frequency  case) 

(HI)  kz  >  k:  kr  is  imaginary  (low  frequency  case) 


Case  I:  Cylindrical  Wave 

When  kz  —  0,  there  is  no  axial  propagation  and  the  impedance  depends  only  on 
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the  acoustic  wavenumber.  In  normalized  form,  it  is: 

pc  H{1]  ( ka ) 

This  may  be  separated  into  a  real  part,  representing  radiation  resistance,  and 
an  imaginary  part,  representing  reactance,  as  follows: 


Re 


J\  (ka)  y0  (ka)  +  Jo  (ka)  Y\  (ka) 
J?(ka)+Y?(ka) 


and 

Y0  (ka)  Yi  (ka)  -  J0  (ka)  Jx  (ka) 

J\  (ka)  +  Y?  (ka) 

These  functions  are  plotted  in  Figure  2-26  as  functions  of  the  non-dimensional 
acoustic  wavenumber,  ka. 

The  radiation  resistance  increases  from  zero  at  ka  =  0  and  asymptotically  ap¬ 
proaches  unity  for  large  values  of  the  argument.  This  may  be  shown  by  considering 


lim  Za  =  —ipc  lim 

ka— >oo  ka^oo 


.7V, 


— ipcex.p(i-~ )  =  pc 
& 


That  is,  the  radiation  resistance  approaches  the  plane  wave  specific  acoustic 
impedance  for  large  values  of  ka. 

The  normalized  radiation  reactance,  Figure  2-26  b,  increases  from  zero  at  ka  = 
0  to  a  maximum  magnitude  of  -0.3971  at  ka  =  0.4617,  subsequently  decreasing 
monotonically  to  zero  at  large  values  of  ka.  The  negative  value  of  the  reactance 
indicates  that  physically  it  is  mass-like,  representing  the  accession  to  inertia  of  the 


external  fluid. 
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Case  II:  High  Frequency 

When  the  axial  wavenumber  kz  is  less  than  the  acoustic  wavenumber  k,  the  radial 


wavenumber  is  real.  The  real  and  imaginary  components  of  the  normalized  spectral 
impedance  are: 


Re 


-i 


(—)  = 

[,  (k*Y] 

2 

J\  (ka)  Yo  (ka)  +  Jo  (ka)  Yi  (ka) 

\pc) 

rwj 

J\  (ka)  +  Yj2  (ka) 

and 


Im 


(—)  = 

[-(«! 

1 

2 

Y0  (ka)  Yi  (ka)  —  J0  (ka)  R  (ka) 

\pc) 

J\  (ka)  +  Y*  (ka) 

Case  III:  Low  Frequency 
For  k  <  kz,th.e  normalized  impedance  is: 


Za=_i 

\(bY  il 

1 

2 

'K0 

(kra) 

pc 

UJ  j 

K, 

^ kra ^ 

Since  K0  (kTaj  and  K\  (kraj  are  both  real,  ~  is  always  imaginary  for  k  <  kz 


For  very  low  frequencies,  i.  e.  k  «  kz,  kr  ~  kz  and  for  kra  »  1, 


KoOh-a) 
K\  (kTa) 


l,and 


lim  Za  =  a  mass  reactance. 

k->  o  kz 

kza»  1 


The  real  and  imaginary  components  of  the  normalized  spectral  impedance  are 
plotted  as  a  function  of  ka  in  Figure  2-27.  These  curves  are  plotted  for  four  values 
of  kza\  |,  7T,  and  2tx. 

The  value  of  kza  divides  the  plots  into  the  low  frequency  (k  <  kz))  and  high 
frequency  (k  >  kz))  regions. 

In  the  low  frequency  region,  the  real  part  of  Za  is  zero  and  there  is  no  radiated 
acoustical  energy.  The  reactance  Im(Za)  increases  with  frequency.  Both  the  real 
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and  imaginary  components  become  infinite  at  the  cut-off  frequency,  where  kr  =  0. 
This  marks  the  transition  between  an  evanescent  radial  wave  component  and  a 
radiating  radial  wave  component. 

In  the  high  frequency  region  the  radiation  resistance  decreases  as  k  increases, 
asymptotically  approaching  pc  when  k  »  kz,  where  kr  k. 

The  radiation  reactance  increases  with  frequency  below  the  cutoff  frequency, 
and  decreases  above  cutoff,  falling  to  zero  at  high  values  of  ka. 
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CHAPTER  3 
Elastic  Media 

3.1  Introduction 

A  concise  review  of  the  historical  development  of  vibration  analysis  is  provided 
by  Soedel  (1981).  Love  (1888)  is  generally  credited  with  having  firmly  established 
the  basis  for  the  analysis  of  vibrations  in  thin  shells.  The  assumptions  made  in 
thin  shell  theory  as  described  by  Junger  and  Feit  (1986)  are:  that  the  thickness 
of  the  shell  is  small  compared  with  any  length  or  surface  dimension  and  with  the 
smallest  radius  of  curvature  of  the  shell;  that  the  displacement  is  small  relative  to 
the  thickness;  that  the  transverse  normal  stress  acting  on  planes  parallel  to  the  shell 
middle  surface  is  negligible;  and  that  fibers  of  the  shell  normal  to  the  middle  surface 
remain  so  after  deformation  and  are  not  subject  to  elongation.  The  problem  is  thus 
reduced  to  the  determination  of  the  deflection  of  the  middle  surface  of  the  shell. 
Although  there  is  no  exact  definition  of  "thinness”,  Krauss  (1967)  suggests  that  a 
maximum  thickness  to  radius  ratio  of  1  to  10  defines  the  limitation  of  thin  shell 
theory.  Herrmann  and  Mirsky  (1956)  made  a  comparison  between  approximate 
shell  theories  with  the  exact  solutions  obtained  using  three  dimensional  elasticity 
theory.  They  concluded  that  when  the  wavelength  of  the  vibration  in  the  shell 
is  greater  than  about  thirty  times  the  thickness  of  the  shell,  simple  membrane 
shell  theory  adequately  describes  the  shell  behavior.  For  wavelengths  greater  than 
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about  eight  times  the  thickness,  shell  bending  theory  yields  results  which  are  in 
good  agreement  with  the  exact  theory. 

Since  the  early  work  of  Love  a  large  number  of  investigators  have  developed 
differing  equations  of  motion  for  thin  cylindrical  shells.  Leissa  (1973)  published 
an  extensive  monograph  in  which  many  of  the  different  thin  shell  theories  are 
compared  in  great  detail.  Leissa’s  monograph  contains  over  500  references  dealing 
with  circular  cylindrical  shells.  Junger  and  Feit  (1986)  emphasize  that  most  of 
the  differences  in  thin  shell  theories  are  related  to  small  differences  in  the  strain- 
displacement  relationships,  which  do  not  significantly  affect  numerical  results. 

The  simplest  theory  for  circular  cylindrical  shells  includes  only  extensional 
(membrane)  effects,  resulting  in  a  system  of  fourth  order  partial  differential  equa¬ 
tions  in  the  axial  and  azimuthal  coordinates.  When  bending  effects  are  included, 
the  resulting  equations  are  of  eighth  order.  The  alternative  formulations  of  the 
’’classical”  shell  theory  result  from  different  approximations  made  in  representing 
the  bending  effects.  When  transverse  shear  and  rotary  inertia  are  included,  a  set 
of  tenth  order  equations  results.  This  case  has  been  investigated  by  Naghdi  and 
Cooper  (1956)  and  by  Lin  and  Morgan  (1956),  using  the  ’’strength  of  materials” 
approach.  These  authors  obtained  the  dispersion  curves  and  amplitude  ratios  for 
the  first  three  axisymmetric  modes  of  an  infinite  cylindrical  shell.  Later  Cooper 
and  Naghdi  (1957)  extended  their  work  to  include  the  non-axisymmetric  modes 
for  the  first  five  axial  modes  with  ten  circumferential  modes. 
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3.2  Thin  Shell  Theory 


3.2.1  Equations  of  Motion  for  Thin  Cylindrical  Shells 

In  the  current  work,  we  shall  employ  the  simplified  equations  of  motion  derived  by 
Donnell  (1933)  as  described  by  Junger  and  Feit  (1986). 

Thus  we  consider  a  thin  cylindrical  shell  in  a  (r,  6 ,  z )  coordinate  system.  The 
thickness  h  of  the  shell  is  small  compared  with  the  mean  radius  a.  The  variables 
ic,u,and  u  represent  the  radial,  circumferential,  and  axial  components  of  dynamic 
displacements  of  the  shell  midsurface,  and  are  assumed  to  be  small  compared  with 
the  thickness  h.  A  constant  (3  is  defined  as 


01  = 


H 2 

12a2 


Terms  proportional  to  (32  in  the  equations  of  motion  are  due  to  the  presence  of 
bending  stresses.  Let  A p  represent  the  differential  pressure  between  the  inside  and 
outside  surfaces  of  the  shell  when  fluid  loading  is  present.  The  Donnell  equations 
then  may  be  written  as: 

Axial: 


d2u  (1  —  v)  d2u  (1  +  v)  d2v  v  dw  1  d2u 
d?  +  2a2  302  +  2a  dzd6  +  ~a!h  ~  = 

Tangential: 

(1  +  v)  d2u  (1  —  v )  d2v  1  d2v  1  dw  1  d2v 

2 a  dzd9  +  2  dz2  a 2  d92  a2  d6  c2  dt 2 


Radial: 


vdu  1  dv  w 

- 1 - 1 - [-  3 

a  dz  a 2  d9  a 2 


2d*w  cftw 
a2—  +  2- 


dz 4  d2z26 


+ 


1  d*w\  ld2w  t  (1-v2) 


a‘ 


d*e  j 


+ 


c2  dt 2 


4-  Ap- 


Eh 
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In  these  equations  E  is  the  Young’s  modulus  of  the  shell  material,  v  is  Poisson’s 
ratio,  and  cv  is  the  low-frequency  phase  velocity  of  compressional  waves  in  an  elastic 
plate,  given  by 

in  which  p  is  the  density  of  the  shell  material. 

3.2.2  Axisymmetric  Case 

In  the  present  work,  we  are  concerned  with  the  axisymmetric  response  of  the 
cylinder.  With  this  restriction,  all  of  the  partial  derivatives  with  respect  to  the 
angular  coordinate  9  vanish  and  the  equations  of  motion  simplify  to: 

Axial: 


d2u  v  dw  1  d2u 
dz 2  a  dz  c2  dt2 

Tangential: 


(3.1) 


(1  —  v)  d2v  1  d2v 
2  dz2  ~  ~\dE 


=  0 


Radial: 


v  du  w  2  2  d^w  1  d2w  (1-v2) 

~7T  +  "2  +  &  a  a T  +  W-5W  +  AP~  El  =  0 
a  dz  a 2  azA  c2  dt2  Eh 


(3.2) 


(3.3) 


We  note  that  in  this  simplified  form,  the  equation  for  tangential  displacement 
v  is  uncoupled  from  the  equations  for  radial  and  axial  motion.  The  latter  two 
equations  are  cross-coupled  via  the  Poisson’s  ratio,  in  the  terms  and  . 
It  is  this  cross  coupling  which  is  responsible  for  most  of  the  complications  in  the 


behavior  of  thin  elastic  shells. 


61 


3.3  Transfer  Matrix  Model  of  a  Thin  Shell 


To  obtain  the  transfer  matrix  for  a  thin  shell,  we  apply  the  following  two- 
dimensional  Fourier  transform  to  the  axisymmetric  thin  shell  equations,  3.1  and 


f(k„w)  =  //  f(z,t)e-‘^+^izd.t 


Application  of  this  transformation  to  3.1  yields: 


_  or  _ 

-kzu  —  iky—w  H — -u  =  0 


a  ct 


and  operation  on  3.3  results  in: 


1  .  o2„2u4  ^ 


—ikz—u  +  w  —  +  /32azk* - -  +  ApT  =  0 


2  C2 


where  u  =  u(kz,u ;),  w  =  w(hz,L>),  and  T  =  ^ 

Solving  for  u  in  3.4,  substitution  into  3.5,  multiplying  through  by  a2,  and 


collecting  terms  results  in 


+  i  +  +Apa2T  =  0 


\  1  ' - ^2 

|  -  n) 


Substitution  of  the  non-dimensional  frequency  Q  =  c oa/cv  and  the  non-dimensional 
axial  wavenumber  —  kza  into  equation  3.6  yields 

If 2?,2  1 

*  (giVj+1+,82y»~tf  +a?“t=o  <3-7) 

The  Fourier  transform  of  the  radial  velocity  is  obtained  from  the  radial  dis¬ 
placement,  Vr(kz,u )  =  —iuw(kz,u>).  Then  the  radial  spectral  impedance  of  the 
thin  shell  is  calculated  as  the  ratio  of  differential  pressure  to  radial  velocity, 

z*(S2-^)  =  d^[(n^f)+1  +  ^-fi2]  <3-8> 
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3.3.1  Two  Port  Matrix  Model 

Now  let  Pi  represent  the  pressure  acting  on  the  inner  wall  of  the  shell,  and  P2 
be  the  pressure  on  the  exterior  surface  of  the  shell.  One  of  the  assumptions  of 
thin  shell  theory  is  that  there  is  no  radial  deformation  of  the  shell,  so  that  the 
radial  velocity  V\  at  the  inner  surface  is  equal  to  the  radial  velocity  V2  at  the  outer 
surface.  The  differential  pressure  across  the  shell  is  A p  =  Pi  —  P2,  and  thus  the 
shell  may  be  represented  by  the  following  two  port  matrix: 


Pi 

1 

P2 

. F1 . 

0  1 

v2 

3.4  Analysis  of  Thin  Shell  in  Wavenumber-Frequency  Space 


(3.9) 


The  study  of  free  wave  propagation  in  thin  cylindrical  shells  has  received  the 
attention  of  many  investigators,  including  the  authors  cited  in  the  introduction  to 
this  chapter. 

Smith  (1955)  used  Kennard’s  equations  of  motion  (Kennard,  1953)  to  compute 
the  phase  velocities  and  relative  displacement  magnitudes  for  all  possible  modes 
of  free  wave  propagation  in  a  thin  shell.  Naghdi  and  Cooper  (1956)  included  the 
effects  of  transverse  shear  deformation  and  rotatory  inertia  in  their  model  of  a  thin 
shell.  Herrmann  and  Mirsky  (1956)  made  a  comparison  of  several  thin  shell  theories 
with  the  exact  results  of  elasticity  theory.  Their  analysis  determined  the  range  of 
wall  thickness  to  axial  wavelength  over  which  the  various  approximate  theories 
are  acceptable  substitutes  for  the  exact  analysis.  For  the  range  of  frequencies  of 
interest  in  this  study,  the  simplified  equations  of  Donnell  (1933)  will  adequately 
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describe  the  behavior  of  free  waves  in  a  thin  cylindrical  shell. 

For  this  analysis  we  make  use  of  equations  3.4  and  3.5,  with  the  pressure  term 
in  3.5  set  equal  to  zero.  Writing  these  equations  in  matrix  form,  we  have: 


(^2  kz'j  ikz  ^ 

u 

0 

-Mzz  (^+/32a2^-f)  _ 

ui 

0 

(3.10) 


In  order  for  free  wave  propagation  to  occur,  the  determinant  of  equation  3.10 
must  equal  zero.  This  is  the  dispersion  equation,  also  known  as  the  characteristic 
equation.  Expanding  the  determinant  of  3.10  and  substituting  the  non-dimensional 
variables  fi  and  Kz  defined  earlier  results  in  the  following  form  for  the  dispersion 
equation: 


fi4  -  a'1  [l  -  K'i  +  0 2K‘]  +  K]  [(l  -  ■  ’)  +  0‘K 


(3.11) 


We  note  that  equation  3.11  is  of  fourth  order  in  Q  and  sixth  order  in  Kz.  The 
only  parameters  which  influence  the  behavior  of  equation  3.11  are  Poisson’s  ratio 
v  and  the  bending  coefficient  /?.  In  this  study  we  shall  use  the  values  v  =  0.37  and 
P  =  0.0103,  which  represent  the  experimental  tube. 

As  discussed  previously,  the  Donnell  theory  is  limited  in  frequency  range.  Using 
the  Herrmann  and  Mirsky  criterion  that  the  minimum  axial  wavelength  is  eight 
times  the  shell  thickness,  we  find  that  the  upper  limit  on  the  non-dimensional  axial 
wavenumber  is  Kz  <  or  Kz  <  7r/15/?\/l2. 

The  roots  of  3.11  express  the  frequency  as  a  function  of  the  wavenumber  for 
specified  values  of  Poisson’s  ratio  and  /?. 
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n 


[1  +  K\  +  /PK*\  ±  ^[1  +  K |  7  A?]2  -  4A?  [(1  -  v2)  +  (PK*\ 

2 


(3.12) 


The  two  roots  of  3.12  define  two  branches  of  the  dispersion  equation  for  the 
n  =  0  mode  of  the  thin  shell.  The  dispersion  curves  for  the  two  branches  are 
shown  in  figure  3-1.  in  which  the  frequencies  are  plotted  as  a  function  of  the 
axial  wavenumber.  The  dimensions  of  the  experimental  tube  were  used  for  these 
calculations. 

Branch  2  is  cut  off  below  the  ’’ring”  frequency,  where  Q  =  Kz,  fr  =  25.1 
kHz,  and  kz  =  1  =  40  /meter.  For  wavenumbers  far  above  the  ’’ring”  frequency, 
Branch  2  is  not  dispersive. 

Branch  1  is  non-dispersive  for  very  small  wavenumbers.  If  (3  —  0,  that  is,  if 
there  is  no  bending  effect,  the  shell  obeys  membrane  theory.  For  this  situation, 
branch  1  will  cut  off  slightly  below  the  ’’ring”  frequency.  With  f3  finite,  the  bending 
terms  serve  to  inhibit  the  cutoff  of  branch  1.  The  wavenumber  at  which  the 
dispersion  curve  of  branch  1  bends  upward  varies  inversely  with  j3,  that  is,  small 
values  of  j3  delay  the  upturn  to  larger  wavenumbers. 

The  phase  velocity,  defined  as  cv  =  ^ ,  is  shown  for  the  two  branches  in  figure 
3-2.  The  phase  velocity  of  branch  2  becomes  infinite  below  the  ’’ring”  frequency, 
and  is  equal  to  the  plate  wavespeed,  cv  =  3939  meters/second  at  high  wavenumbers. 
The  phase  velocity  of  branch  1  approaches  the  longitudinal  bar  velocity,  = 
cp\/l  —  v2  =  3659  meters/second  for  very  low  wavenumbers  (Graff,  1975).  As  the 
wavenumber  increases,  the  phase  velocity  of  branch  1  decreases  to  a  minimum  and 
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then  begins  to  rise  again.  The  minimum  phase  velocity  for  branch  2  is  given  by 
Cipnun  =  (Herrmann  and  Mirsky,  1956),  and  for  the  experimental  tube, 

c^min  =  544  meters/second.  Lin  and  Morgan  (1956)  developed  an  approximate 
equation  for  the  frequency  of  minimum  phase  velocity.  Their  result  may  be  put 
into  the  form  /min  =  fr^2  (1  -  v2).  where  fr  represents  the  ring  frequency.  For 
the  experimental  tube,  fm-m  =  33  kHz. 

The  group  velocity  is  defined  as  cg  =  ^£12  and  is  plotted  for  the  two  branches 
in  figure  3-3.  The  group  velocity  for  branch  2  increases  from  zero  as  the  branch 
cuts  on  at  the  ring  frequency,  then  levels  off  at  the  plate  wavespeed  cp  at  high 
frequencies.  The  group  velocity  for  branch  1  is  equal  to  the  bar  velocity  c*,  at  low 
wavenumbers  and  drops  rapidly  as  the  ”  ring”  frequency  is  approached.  A  minimum 
of  46  meters/second  occurs  at  kz  =  121  /meter,  corresponding  to  a  frequency  of 
about  23.3  kHz.  From  this  minimum,  the  group  phase  velocity  gradually  increases 
with  wavenumber  and  frequency. 

Figure  3-4  presents  the  ratio  of  radial  to  axial  displacement  for  each  branch  as  a 
function  of  the  axial  wavenumber.  This  figure  indicates  that  at  low  wavenumbers, 
the  displacement  for  branch  2  is  primarily  radial,  and  that  at  high  wavenumbers, 
the  displacement  for  branch  1  is  primarily  radial. 

3.5  Modes  in  Harmonically  Driven  Finite  Length  Tube 

In  order  to  compare  the  experimental  results  with  those  predicted  from  theory, 
it  is  essential  to  accurately  know  the  value  of  wave  speed  in  the  experimental  tube. 
This  problem  is  discussed  in  detail  in  Chapter  5.  The  most  precise  method  was 
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based  on  measurements  of  the  frequency  of  longitudinal  modes  when  the  tube 
was  driven  with  a  continuous  sinusoidal  excitation.  In  this  section  we  develop  the 
theory  used  for  the  wave  speed  determination. 

Consider,  then,  a  thin  cylindrical  shell  of  length  L,  thickness  h,  and  mean  radius 
a,  which  is  freely  suspended  and  driven  radially  by  an  axisymmetric  ’’ring”  source 
at  one  end.  Since  the  excitation  is  axisymmetric,  the  tangential  displacement  v 
is  uncoupled  from  the  radial  displacement  w  and  the  axial  displacement  u.  We 
shall  use  the  axial  and  radial  equations  for  axisymmetric  motion,  3.1  and  3.2.  For 
frequencies  at  which  the  axial  wavelength  \z  is  sufficiently  large,  that  is,  Xz  > 
30 h,  the  term  involving  62  may  be  neglected.  This  condition  is  satisfied  for  the 
experimental  tube  for  frequencies  below  the  ring  frequency,  about  25  kilohertz. 
The  resulting  equations  are  identified  as  those  for  a  thin  membrane,  viz: 


Axial: 


d2u  v  dw 
dz2  ^  a  dz 


1  d2u 
c2  dt2 


(3.13) 


Radial: 


v  du  w  1  d2w 
a  dz  +  a2  c2  dt2 


(3.14) 


The  ’’ring”  displacement  input  at  z  =  0  provides  one  boundary  condition,  viz 


w( 0,  t )  =  Wq  •  exp(-iut) 


(3.15) 


The  second  boundary  condition  is  found  from  the  requirement  that  the  axial 
strain,  must  be  zero  at  the  end  of  the  tube,  where  z  =  L. 

Since  a  harmonic  input  is  assumed,  equations  3.13  and  3.14  are  reduced  to 
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ordinary  differential  equations.  The  axial  displacement  is  described  by 


d2u  v  dw  (  u)  \ 

d?  +  ad^+[Zl 


(3.16) 


and  the  radial  displacement  (with  the  bending  term  deleted)  is  given  by 

v  du  1  (  uj  \ 2  _  /o  i  -7\ 

~—+w  -r  -  -  =0  (3.17) 

adz  az  \cP  J 

Re-arranging  equation  3.17  and  differentiating  with  respect  to  z  yields 

dw  _  va  d?u  /Q  n 

dz  ~  (Q2  -  1)  dz>  {  '  } 

where  the  non-dimensional  frequency  Q.  =  ^  has  been  introduced.  Substitu¬ 
tion  of  3.18  into  3.16  results  in 

d2u 


dz 2 


+  rfu  =  0 


(3.19) 


where 


2  n2(o?-i) 

7  “  a2  [fi2  -  (1  -  i/2)] 


(3.20) 


B, 


Now  we  assume  a  standing  wave  solution  to  equation  3.19, 

u(z)  =  A  cos  7 z  +  B  sin  7 z  (3.21) 

The  boundary  condition  w  =  Wo*exp(—iujt)  at  z  =  0  determines  the  coefficient 


and  the  boimdary  condition  ^  =  0  at  z  =  L  results  in 


A  —  B  cot  7  L 
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Using  these  results,  the  axial  and  radial  displacements  are  given  by 


COS 


7  (L-z) 


sin  7  L 


(3.22) 


and 


w(z)  =  W0 


sin  7  (L  —  z) 


sin  7  L 


(3.23) 


Modal  resonances  will  occur  for  values  of  7  for  which  sinyL  =  0,  or  7  = 
^ ,  m  =  1,  2,  3  ... .  Note  that  the  axial  displacement  is  zero  at  Cl  =  1,  the  ring 
frequency,  and  also  at  Q,  =  y/1  —  z/2,  which  for  brass  (z/  =  .37)  is  at  0.93  times  the 
ring  frequency. 

The  non-dimensional  modal  frequencies  may  be  found  from  the  real  roots  of 


U4-U2 


/  TO7TG\2 

1+(— ) 


/ m7ra\  2  /  9\ 

+  (— )  (i-"2)"0 


or 


(3.24) 


The  expected  value  of  modal  frequencies  is  then  calculated  for  an  assumed 
value  of  extensional  wave  speed  cp  for  the  non-dimensional  frequencies  in  equation 
3.24  using 


c  Q 

/<  _ 

rn  —  ~Z 

Z7 TCI 


(3.25) 


These  are  the  relationships  needed  for  the  wave  speed  determination  described 


in  Chapter  5. 
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3.6  Membrane  Transition  Frequencies 

It  is  convenient  at  this  point  to  develop  an  aspect  of  membrane  theory  which 
vail  be  of  use  in  the  next  topic.  Thus  we  consider  a  harmonically  driven,  semi¬ 
infinite  shell  which  meets  the  criteria  for  analysis  as  a  membrane,  as  discussed 
previously.  As  before,  the  shell  is  radially  driven  at  z  —  0  with  a  harmonic  input 
described  by  equation  3.15.  However,  for  the  present  analysis,  we  shall  assume  a 
travelling  wave  solution  for  the  axial  displacement  u  of  the  form: 

u(z,t )  =  Aexp[— i(qz  +  cut)]  +  Bexp[i(7.z  —  cut)]  (3.26) 

Since  we  have  assumed  a  semi-infinite  tube  extending  in  the  positive  z  direc¬ 
tion,  there  can  be  no  wave  travelling  in  the  negative  z  direction,  and  therefore 
A  —  0.  The  radial  displacement  w  is  given  by 

w(z,  t)  =  —  -----hB  exp[z(7Z  -  wt)\  (3.27) 

Substituting  3.15  into  3.27  and  evaluating  at  z  =  0  determines  B, 

B  w.ffi2  - 1) 

ivaj 

Thus  the  equations  of  motion  for  this  case  are: 

u(z,  t)  =  - —  exp[?(72  —  cut)]  (3.28) 

iva'y 

and 


w(z ,  t )  =  Wq  exp[i(72  —  cut)] 
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Replacing  7  in  the  denominator  of  equation  3.28  with  the  expression  in  equa¬ 
tion  3.20  and  simplifying  results  in 


,  w'ov'tn2  - 1)  -  (1  -  «2)) 

u(z,t)  = - 1 - — ^ - exp [*(72  -  ut) ] 


(3.29) 


Equation  3.29  implies  that  the  horizontal  displacement  u  will  vanish  for  two 
values  of  the  non-dimensional  frequency,  Cl.  These  are  Cl  =  1,  which  is  the  well 


known  ring  frequency,  and  Cl  —  \/l  —  "f2,  which  we  shall  refer  to  as  the  ’’lower  tran¬ 
sition  frequency.”  For  free  waves  propagating  in  an  infinite  membrane  shell,  this 
lower  transition  frequency  corresponds  to  the  cutoff  of  the  flexural  branch  (iden¬ 
tified  as  ’’Branch  1”  earlier).  The  ring  frequency,  or  upper  transition  frequency, 
corresponds  to  the  onset  of  propagation  in  the  longitudinal  mode. 

Now  consider  the  wave  number  7  for  the  present  case.  From  equation  3.20, 

Cl 


7  — 


Cl2-  1 


(3.30) 


a  ^  Cl2  -  (1  -  v2) 

If  there  were  no  Poisson  coupling,  that  is,  v  =  0,  then  7  =  Cl /a  =  u/cp,  which  is 
identified  as  an  ordinary  wave  number  for  a  non-dispersive  wave.  For  values  of  SI  < 
VT  —  v2  and  for  Cl  >  1,  the  expression  under  the  radical  in  equation  3.30  is  positive, 


and  7  is  real.  For  yf\--v2  <Cl  <  1,  the  expression  is  negative,  and  7  is  imaginary. 
Thus,  7  ’’jumps”  from  a  real  number  to  an  imaginary  number  at  the  lower  transition 
frequency,  and  back  to  a  real  number  at  the  upper  transition  (or  ’’ring”)  frequency. 
For  the  brass  tube  used  in  this  study,  the  lower  transition  frequency  is  23.338 
kilohertz,  and  the  upper  transition  frequency  is  25.121  kilohertz.  Also  note  that  7 


has  a  pole  at  Cl  —  \ZT--~v2,  and  a  zero  at  Cl  =  1. 
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3.7  Transient  Response  of  Thin  Shell 

3.7.1  Previous  Work  on  Transient  Response 

In  1948,  R.  M.  Davies  published  an  82  page  paper  titled  ”A  Critical  Study  of 
the  Hopkinson  Pressure  Bar.”  (The  Hopkinson  pressure  bar  is  an  apparataus 
which  was  used  to  study  the  propagation  of  stress  pulses;  a  concise  description  is 
given  by  Kolsky  (1963)).  Davies  used  the  exact  equations  of  motion  of  a  circular 
elastic  cylinder  developed  by  Pochhammer  (1876)  and  Chree  (1889),  as  presented 
by  Love  (1944)  to  compute  the  dispersion  curves  (phase  and  group  velocity  versus 
wavenumber)  for  the  first  three  modes  of  vibration  of  an  infinite  length  bar.  He 
used  the  method  of  stationary  phase  to  predict  the  displacements  and  stresss  in 
a  free  bar  excited  by  an  impulse  applied  one  end.  The  group  velocity  curves  of 
the  first  two  modes  exhibit  a  ”  valley”  which  is  also  a  characteristic  of  the  flexural 
mode  of  an  axysymmetric  thin  shell. 

Using  a  very  ingenious  physical  argument,  Skalak  (1957)  constructed  an  in¬ 
tegral  solution  for  the  propagation  of  the  pulse  resulting  from  the  longitudinal 
impact  of  two  semi-infinite  rods.  The  integral  was  evaluated  via  the  asymptotic 
saddle-point  method  for  large  values  of  time.  The  resulting  pulse  is  described 
by  Airy  integrals  and  is  characterized  by  a  smooth  rise  from  zero  to  a  maximum 
value,  followed  by  a  decaying  oscillation  about  the  level  which  is  predicted  from 
the  elementary  (non-dispersive)  theory  using  the  one  dimensional  wave  equation. 

Berkowitz  (1963)  analyzed  the  response  of  a  semi-infinite  cylindrical  shell  mov¬ 
ing  at  constant  velocity  which  impacts  a  rigid  surface.  The  impact  provides  an 
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axisymmetric  step  function  velocity  input  in  the  longitudinal  direction.  He  for¬ 
mulated  the  problem  in  terms  of  the  Laplace  transformed  (in  time)  axysimmetric 
membrane  equations  for  the  radial  and  longitudinal  displacements,  and  the  axial 
and  circumferential  stresses. 

The  resulting  inversion  integrals  cannot  be  evaluated  in  closed  form,  so  the 
method  of  stationary  phase  was  used  to  obtain  approximate  asymptotic  solutions 
for  large  values  of  non-dimensional  time.  The  solution  is  characterized  by  an 
initial  disturbance  which  propagates  at  the  plate  wavespeed,  cp,  followed  by  a 
steep  wavefront  travelling  at  the  bar  velocity,  c0.  This  is  followed  by  a  damped 
oscillation  which  is  described  by  the  Airy  integral. 

An  excellent  discussion  of  these  and  other  investigators  of  transient  phenomena 
is  presented  by  Graff  (1975). 

3.7.2  Synopsis  of  Approach 

The  objective  of  the  present  study  is  to  compute  the  predicted  response  of  the 
experimental  cylindrical  shell  to  an  axisymmetric  radial  impulse  input  applied  at 
the  origin  of  the  axial  coordinate,  z  =  0.  A  detailed  description  of  the  experiment 
is  presented  in  Chapter  5. 

The  axisymmetric  equations  of  motion  for  the  radial  and  axial  displacements 
(Equations  3.1  and  3.3)  are  Fourier  transformed  in  time  and  Laplace  transformed 
in  the  axial  coordinate  z.  The  resulting  transformed  equations  are  algebraic  in  the 
transform  variables  kp  and  s  where  kp  —  u>/cp  is  a  frequency  dependent  wavenum¬ 
ber,  and  5  is  a  complex  spatial  wavenumber,  s  =  a  +  ikz ,  where  kz  is  the  real 
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valued  axial  wavenumber.  An  axisymmetric  radial  delta  function  input  is  applied 
to  the  system  at  the  origin,  z  =  0.  The  inverse  Laplace  transform  of  the  result  is 
then  computed  by  contour  integration,  which  yields  two  equations  in  ( kp ,  z)  space. 
These  equations  are  essentially  transfer  functions  which  describe  the  radial  and 
longitudinal  displacements  at  axial  location  z  as  a  function  of  frequency  lj,  that  is, 
u(u,  z )  and  w  (u,  z ).  When  these  functions  are  multiplied  by  the  Fourier  transform 
of  an  input  at  z  =  0,  and  the  inverse  Fourier  transform  is  applied  to  the  product, 
the  result  is  the  response  as  a  function  of  time  at  the  axial  position  specified;  that 
is,  u(t,z0 )  and  w  ( t,z0 )  where  z0  is  a  particular  observation  point. 

3.7.3  Transformed  Equations  of  Motion 

A  double  transform  comprising  a  Fourier  transform  on  the  time  variable,  and  a 
Laplace  transform  on  the  axial  variable  z  defined  as: 

^  roc  roc 

f{kp,s)=  /  eikpt  /  f(t,  z)e~S2dz  dt 

J  oo  J0 

is  applied  to  the  axysimmetric  equations  of  radial  and  longitudinal  motion,  3.1 
and  3.3  .  The  resulting  transformed  equations  of  motion  are: 

Axial: 


S(s2  +  fcJ)+™(s-)  =  0  (3.31) 

Cl 

Radial: 


u(s-)  +  w(\  +  s4/32 a2  -  k2p)  =  P{u,  s ) 
a  a1 


(3.32) 
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where  P(u,  s )  is  the  transform  of  the  input, 

^  roc  roc 

P(w,s)=  P(t,z)e-szdzdt 

Joe  Jo 

The  mechanical  input  for  the  experiment  is  provided  by  a  piezoelectric  disk  which 
is  bonded  to  the  end  of  the  brass  tube.  For  transient  measurements,  a  short 
rectangular  pulse  of  voltage  is  applied  to  the  disk.  This  results  in  a  radial  expansion 
of  the  disk,  which  in  turn  applies  an  axisymmetric  radial  displacement  to  the  tube. 
We  may  approximate  this  input  as  P0S(t)  S(z),  that  is,  a  delta  function  of  time 
and  axial  position,  in  which  P0  has  the  dimensions  of  inverse  length.  Applying  the 
double  transform  above  to  this  input  yields: 

P(uj,  s )  =  P0 

Substitution  of  this  input  into  equation  3.32  and  solving  for  u  and  w  results  in 

Q0  1p 

u(kr,s)  =  --^  (3.33) 

and 

w(kp,s)=(-s2  +  **)P°  (3.34) 

In  equations  3.33  and  3.34,  D  represents  the  dispersion  equation,  which  is  a 
polynomial  of  sixth  order  in  s  and  fourth  order  in  kp,  as  indicated  below: 


D  =  s\(3af  +  s\(3akpf  +  s2 


(1  -  u2) 


-K 


az 


a 


(3.35) 
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3.7.4  Response  in  Spatial  and  Frequency  Domain 

To  obtain  the  radial  and  longitudinal  displacements  as  functions  of  axial  position 
2:  and  frequency  (or  wavenumber,  kp),  we  apply  the  complex  inversion  integral 
(Spiegel,  1971)  to  equations  3.33  and  3.34, 

1  rc+ioo 

9(z)=t~-  e“G(s)ds,  z>0 

£7TZ  Jc— zoo 

This  is  accomplished  by  contour  integration,  using  the  Cauchy  residue  theorem 
(Hildebrand,  1962),  that  is, 

<f  eszG(s )  ds  =  2m  £  R(ak)  (3.36) 

k=l 

where  R(ak)  is  the  A:th  residue.  The  residues  are  evaluated  at  the  singularities 
of  G(s). 

We  note  that  u(kp,s)  and  w(kp,  s)  share  the  same  singularities,  namely,  the 
roots  of  equation  3.35. 


Location  of  Roots 


Consider  the  expression  for  D,  equation  3.35.  Although  D  is  of  sixth  order  in  s, 
it  is  a  cubic  in  s2.  Now  let  q  =  s2;  so  that  D  may  be  written  as 


D  =  qs{/3( x)2  +  q2(/3akP)2  +  q 


ft-*2) 

2  KV 


a * 


(3.37) 


Equation  3.37  will  have  three  roots.  The  roots  will  either  be  all  real,  or  one  real 
and  one  complex  conjugate  pair,  depending  on  the  discriminant  of  3.37.  (Spiegel, 
1968)  Denoting  the  roots  of  3.37  as  Qn,  n  =  1,  2,  3,  the  six  roots  of  equation 
3.35  are  given  by: 

Sj  =  ±yfo'n,  j  =  6  (3.38) 
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The  values  of  S  may  be  real,  imaginary,  or  complex,  depending  on  the  values 
of  Q.  All  of  these  roots  will  be  functions  of  frequency  due  to  the  variable  kp  =  u/cp 
in  equations  3.35  and  3.37.  Although  exact  literal  expressions  for  the  values  of  the 
roots  S  may  be  written,  the  results  are  so  cumbersome  algebraically  that  little 
insight  into  the  behavior  of  the  roots  as  a  function  of  frequency  is  achieved.  Thus 
in  this  study  we  pursue  a  numerical  evaluation  using  the  following  values  of  the 
parameters  appropriate  for  the  experimental  tube: 

Mean  radius,  meters:  a  —  .0249555 

Thickness,  meters:  h  =  .000889 

Poisson’s  ratio:  v  =  0.37 

Plate  wave  speed,  meters/second:  cp  =  3939 

The  value  of  cp  was  determined  experimentally,  as  discussed  in  Chapter  5.  The 
frequency  range  of  interest  in  this  study  is  from  0  to  100  kilohertz,  about  four  times 
the  ’’ring”  frequency,  fr  =  cp/27ra,  which  is  25.12  kHz  for  the  experimental  tube. 


Variation  of  roots  Qn  with  Frequency: 

We  consider  the  nature  of  the  roots  Qn  of  equation  3.37  in  three  frequency 
regimes:  (1)  far  below  the  lower  transition  frequency;  (2)  far  above  the  upper 
transition  frequency;  and  (3)  in  the  range  between  the  transition  frequencies.  The 
real  and  imaginary  components  of  the  three  roots  are  shown  in  figure  3-5  for  the 
frequency  range  from  10  kilohertz  to  40  kilohertz. 
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(1)  Below  lower  transition  frequency:  Below  the  lower  transition  frequency, 
roots  Q\  and  Q%  are  a  complex  conjugate  pair.  Root  Qz  is  negative  and  real. 

(2)  Above  the  upper  transition  frequency:  Far  above  the  ’’ring”  frequency, 
Q\  and  Qz  are  negative  and  real,  and  Q%  is  positive  and  real. 

(3)  Between  the  transition  frequencies:  In  this  range,  the  behavior  of  the 
roots  becomes  very  complicated.  Figure  3-6  presents  the  frequency  dependence 
of  the  roots  in  the  region  between  the  transition  frequencies.  Root  Q i,  which 
is  complex  with  a  positive  imaginary  component  at  low  frequencies,  experiences 
an  increase  in  the  real  component  slightly  below  the  first  transition  frequency. 
The  imaginary  component,  which  has  been  decreasing  with  increasing  frequency, 
suddenly  drops  to  zero  at  the  lower  transition  frequency.  The  real  part  of  Q\ 
jumps  from  a  positive  to  a  negative  value,  and  continues  to  become  more  negative 
at  higher  frequencies. Root  is  complex  with  a  negative  imaginary  component 
at  low  frequencies.  The  imaginary  component  jumps  to  a  positive  value  at  the 
lower  transition  frequency,  then  decreases  to  zero.  The  positive  component  of  Qz 
continually  increases  with  frequency. 

Root  Qz  is  negative  and  real  at  low  frequencies,  becoming  more  negative  as 
the  lower  transition  frequency  is  approached,  at  which  point  it  becomes  complex 
for  a  small  range  of  frequencies.  The  real  component  jumps  to  a  positive  value 
at  the  lower  transition  frequency  and  remains  positive  up  to  the  upper  transition 
frequency.  We  shall  find  shortly  that  the  bending  term  /3  has  an  important  effect 
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on  these  roots  in  the  transition  region. 

Figure  3-7  shows  the  behavior  of  the  magnitude  of  all  three  roots  in  the  region 
between  the  transition  frequencies.  Below  the  lower  transition  frequency,  roots 
Q i  and  Q2  are  equal  in  magnitude  and  are  decreasing.  Root  Qz  (dashed  line)  is 
increasing.  All  three  roots  meet  at  the  lower  transition  frequency,  23.338  kHz.  At 
this  frequency,  root  Q 1  becomes  real  and  negative,  increasing  in  magnitude  with 
frequency.  Roots  Q2  and  Qz  now  form  a  complex  conjugate  pair  until  a  second 
branch  point  is  reached  at  23.606  kHz.  At  the  second  branch  point,  Q2  becomes  real 
and  positive,  increasing  in  magnitude  with  frequency.  Root  Qz  becomes  real  and 
positive,  decreasing  rapidly  with  frequency  until  the  ’’ring3’,  or  upper  transition, 
frequency,  where  it  vanishes.  Above  the  ’’ring”  frequency,  Qz  is  real  and  negative, 
increasing  in  magnitude  with  frequency.  The  branch  points  just  described  are  illus¬ 
trated  in  figure  3-8.  where  the  path  of  the  Q  roots  is  plotted  in  the  neighborhood 
of  the  lower  transition  frequency. 

Roots  of  Dispersion  Equation 

The  six  roots  of  equation  3.37  in  the  complex  wavenumber  (s)  plane  are  found 
from  the  values  of  Qn.  as  indicated  in  equation  3.38.  Negative  real  values  of  the 
roots  Qn  give  rise  to  imaginary  values  of  Sj.  The  six  roots  in  s  occur  in  conjugate 
pairs,  from  taking  both  the  positive  and  negative  values  of  %/Q^,.  The  relationship 
between  the  Q  roots  and  the  S  roots  is  shown  qualitatively  in  figure  3-9.  Here 
the  relative  positions  of  the  roots  are  shown  in  the  complex  Q  and  S  planes.  The 
arrows  indicate  the  path  followed  by  the  roots  as  the  frequency  is  increased.  In 
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figures  3-9  (a  &  b)  the  paths  of  the  roots  are  shown  for  frequencies  below  the  lower 
transition  frequency.  Figures  3-9  (c  &  d)  present  the  relationships  for  frequencies 
above  the  upper  transition  frequency. 

Variation  of  Roots  Sn  With  Frequency:  For  brevity,  we  shall  focus  on  a  de¬ 
scription  of  the  behavior  of  the  S  roots  which  arise  from  the  positive  square  root 
of  the  Q  roots.  A  similar  argument  is  easily  constructed  for  the  image  roots. 

Below  the  lower  transition  frequency,  roots  Si  and  S2  are  a  complex  conju¬ 
gate  pair  with  positive  real  part.  The  magnitude  of  these  roots  decreases  as  the 
frequency  increases.  Root  S3  derives  from  root  Q3,  which  is  negative  and  real  in 
this  frequency  regime.  Therefore,  S3  is  positive  and  imaginary  below  the  lower 
transition  frequency.  Root  S3  increases  with  frequency. 

Above  the  upper  transition  frequency,  roots  Si  and  S3  are  both  positive  and 
imaginary.  These  roots  increase  with  frequency.  Root  S2  is  positive,  real,  and  also 
increases  with  frequency. 

In  the  region  between  the  lower  and  upper  transition  frequencies,  the  mag¬ 
nitudes  of  the  S  roots  vary  in  a  fashion  similar  to  the  behavior  of  the  0  roots 
described  above.  This  is  shown  in  figure  3-10.  The  changes  in  the  root  S3  are  of 
special  interest  and  will  be  analyzed  further. 

Figure  3-11  (a)  is  a  plot  of  the  real  and  imaginary  components  of  S3  in  the 
transition  zone.  Note  that  this  root  acquires  a  real  component  in  this  frequency 
range.  This  is  the  only  part  of  the  entire  spectrum  for  which  root  S3  has  a  real 
part.  The  imaginary  part  of  S3  is  positive  below  the  lower  transition  frequency 
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(identified  by  the  symbol  ”fl”  on  figure  3-11  (a)),  where  it  ” jumps”  to  a  negative 
value.  It  remains  negative  until  the  second  branch  point  (symbol  ”f2”)  at  23.606 
kHz,  where  it  vanishes.  The  imaginary  component  continues  to  be  zero  until 
the  ’’ring”  frequency  (symbol  ”f3”).  At  this  frequency,  the  real  component  of  S3 
vanishes,  and  the  imaginary  component  returns,  becoming  positive  and  increasing 
with  frequency. 

The  rather  unusual  behavior  of  the  S  roots  in  the  transition  zone  is  related  to 
the  bending  term  (3  in  the  governing  equation  for  the  radial  diplacement,  equation 
3.3.  This  will  be  demonstrated  by  plotting  the  components  of  S3  in  the  transition 
zone  for  several  values  of  f3. 

Figures  3-11  (b),  3-11  (c),  and  3-11  (d)  present  the  components  of  root  S3  for 
values  of  (3  =  .1,  .001,  and  .0001,  respectively.  We  see  that  as  (3  decreases,  the 

second  branch  frequency  moves  towards  the  lower  transition  frequency.  In  addition, 
the  magnitude  of  S3  increases  substantially  at  the  lower  transition  frequency.  In 
the  limit  as  (3  — >  0,the  real  part  of  S3  and  the  second  branch  frequency  both  vanish. 
This  limit  represents  the  transition  to  shell  membrane  theory,  for  which  there  will 
only  be  two  roots. 

Specification  of  Contour 

In  order  to  apply  the  Cauchy  residue  method  indicated  in  equation  3.36,  the 
contour  of  integration  must  be  specified.  Since  we  are  interested  in  wave  propaga¬ 
tion  in  the  positive  z  direction,  we  choose  a  Bromwich  contour  which  includes  the 
imaginary  axis,  and  a  semi-circle  of  infinite  radius  in  the  left  half  of  the  complex 
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s  plane.  This  contour  is  shown  in  figure  3-12  for  the  pole  locations  in  the  low  fre¬ 
quency  regime  (3-12  (a))  and  for  the  high  frequencies  (3-12  (b)).  We  include  the 
poles  which  derive  from  the  negative  roots  of  Qn,  that  is,  S4,  S5,  and  S6.  Note  that 
the  poles  which  are  excluded  represent  waves  travelling  in  the  negative  2  direction. 

Evaluation  of  Residues 


We  now  express  equations  3.33  and  3.34  in  terms  of  the  poles  which  are  to  be 
included,  viz. 

svP0 


u(kp,s) 


and 


w 


(kp,  s ) 


'a(s-SA)  ( s-S5 )  (s  -  S6) 


(s2  +  kfj  Pp 


(3.39) 


(3.40) 


(s  -  S4)  (s  -  S5 )  (s  -  S6 ) 

Each  of  equations  3.39  and  3.40  will  have  three  residues,  which  we  shall  desig¬ 
nate  as  UR1,UR2,UR3  and  WRl,WR2,WR3,  respectively.  Now  if  the  function 
G(s )  in  equation  3.36  has  a  simple  pole  at  s  =  a*,  for  instance,  the  fcth  residue  is 
found  as: 


R(k)  =  lim  [(s  —  afc)  G(s )  es 


Applying  equation  3.41  to  3.39  and  3.40,  we  have: 


(3.41) 


UR1  = 


vP0S4es*z 


a(S4  -  S5){S4  -  S6) 


(3.42) 


UR2  = 


v  Pq  S5  es$z 


a(S5-Si)(S5-Se) 


(3.43) 
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and 


UR3  = 


v  P0  S6  eSeZ 


WR1  - 


WR2 


WR3  = 


a{S6 

-  St)(So  - 

-St) 

Po 

(s? + *!) 

eS4z 

'  (S* 

-  Soils,  - 

-S6) 

Po 

{Si +  4) 

gS5z 

'  (5s 

—  Si)(S$  - 

-S6) 

Po 

(si  +  %) 

es« 

(3.44) 


(3.45) 


(3.46) 


(3.47) 


(S,  -  SiKSe  -  St) 

Equations  3.42  to  3.47  are  essentially  transfer  functions  which  describe  the 
frequency  dependence  of  the  three  components  of  longitudinal  displacement  u  and 
the  three  components  of  radial  displacement  w  as  functions  of  the  axial  position, 
z.  The  net  displacements  are  obtained  by  summing  the  components,  viz. 


w(kp,  z )  =  WR1  +  WR2  +  WR3  (3.48) 


u(kp,  z)  =  UR1  +  UR2  +  UR3 


(3.49) 


3.8  Discussion  of  w(f,  z}  and  u(  f,  z)  Transfer  Functions 


We  now  consider  the  transfer  functions  for  the  radial  and  longitudinal  displace¬ 
ments  in  the  spatial  (z)  and  frequency  (/  =  cpkp/2?c)  domains.  Figures  3-13  to 
3-16  present  three  dimensional  surface  plots  of  the  magnitude  (absolute  value)  of 
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the  displacement  functions  w(f,z)  and  u(f7z)  for  frequencies  ranging  from  0  to 
50  kilohertz,  and  two  ranges  of  axial  distance,  0  to  5  cm,  and  1  to  50  cm.  Plots 
with  a  logarithmic  scale  for  the  functions  are  also  included.  The  logarithmic  plots 
enhance  low  amplitude  characteristics  which  are  not  obvious  on  a  linear  scale.  A 
small  offset  (.001)  is  added  to  the  magnitude  of  the  functions  to  provide  a  ’’floor” 
for  the  logarithmic  plots. 

The  far-field  topography  of  w(f:z)  is  emphasized  in  figure  3-13  (a),  which 
exhibits  a  plateau  above  the  ’’ring”  frequency,  a  steep  cliff  in  the  transition  region, 
and  a  rather  flat  plain  below  the  lower  transition  frequency.  There  is  a  ripple  near 
the  edge  of  the  plateau  which  indicates  a  variation  in  the  magnitude  of  the  spectral 
components  with  axial  distance.  The  logarithmic  plot  of  this  function  (figure  3-13 
(b))  shows  a  smooth  rise  from  the  plain  to  the  plateau.  The  ripples  at  the  edge  of 
the  plateau  are  suppressed  in  this  presentation.  The  ’’wall”  at  z  =  0  cm  is  due  to 
the  near  field. 

The  predominant  feature  in  the  far  field  of  the  longitudinal  displacement 
u(f,  z),  figure  3-14  (a),  is  a  jagged  ridge  which  is  parallel  to  the  z  (distance)  axis. 
This  ridge  is  located  in  the  transition  region  in  frequency  space..  The  function 
falls  off  smoothly  at  lower  frequencies.  The  logarithmic  plot,  figure  3-14  (b),  shows 
a  number  of  smaller  ridges  which  radiate  away  from  the  primary  ridge  above  the 
transition  region.  These  minor  ridges  follow  a  curved  trajectory  in  (f,z)  space. 

Figures  3-15  and  3-16  focus  on  the  near  field.  The  ’’wall”  parallel  to  the 
frequency  axis  at  z  =  0  in  figure  3-15  (a)  reflects  the  fact  that  the  temporal  input 
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to  the  model  is  a  delta  function  in  the  time  domain,  which  has  a  uniform  spectral 
density  in  the  frequency  domain.  The  presence  of  evanescent  near  field  frequency 
components  is  evident  near  the  origin  of  figure  3-15  (a).  The  near  field  stucture  of 
w(f,z)  is  more  visible  in  the  logarithmic  plot,  figure  3-15  (b). 

The  near  field  of  u(f,z)  is  shown  in  figures  3-16  (a)  and  3-16  (b).  The  am¬ 
plitude  of  u  (/,  z)  is  zero  for  all  frequencies  at  u  =  0,  since  the  axial  strain,  is 
zero  at  this  location.  The  formation  of  the  ridge  in  the  transition  zone  as  the  axial 
distance  increases  is  clearly  shown  in  this  figure.  The  minor  ridges  which  radiate 
in  (/,  z)  space  are  visible  in  the  logarithmic  plot,  figure  3-16  (b). 

The  detailed  structure  of  w(f,z)  and  u(f,z )  is  studied  by  taking  ’’slices”  of 
these  functions  at  selected  valaues  of  /  and  z.  Figure  3-17  (a)  is  a  plot  of  the 
magnitude  of  w  (/,  z)  as  a  function  of  frequency  at  four  axial  positions  in  the  near 
field.  At  z  =  .01  cm,  the  spectrum  is  nearly  flat.  This  is  the  upper  edge  of  the 
’’wall”  seen  in  figure  3-15  (a).  As  the  distance  increases,  the  lower  frequencies 
diminish,  whereas  the  high  frequencies  persist.  The  transition  zone  ’’ridge”  is 
evident  in  this  figure,  at  z  =  .11  and  z  =  .21  cm.  Figure  3-17  (b)  shows  the 
frequency  spectrum  at  z  =  .91,  .41,  and  2.01  cm,  indicating  the  formation  of  the 
’’plain”  and  ’’plateau”  regions  as  the  far  field  is  approached.  In  figure  3-17  (c),  the 
frequency  spectrum  is  plotted  for  .2,  .4,  .6,  and  50  cm.  The  far  field  structure  is 
indicated  by  the  spectrum  at  z  =  50  cm. 

The  development  of  the  topography  of  u(f,  z)  is  illustrated  in  figure  3-18. 
Figures  3-18  (a)  and  3-18  (b)  show  the  growth  of  the  ridge  in  the  transition  zone 
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of  u(f,z).  The  spectrum  at  z  =  50  cm  in  figure  3-18  (c)  indicates  the  formation 
of  the  minor  ridges  in  frequency  space. 

The  magnitudes  of  w(f,  z)  and  u(f,  z)  as  a  function  of  axial  distance  is  plotted 
in  figures  3-19  to  3-22  for  several  frequencies  of  interest.  The  near  field  plots, 
figures  3-19  and  3-20  show  the  decay  of  the  low  frequency  components  at  10  and 
22  kilohertz  as  the  far  field  is  approached. 

Figures  3-21  and  3-22  emphasize  the  far  field  characteristics  of  w(f,z )  and 
u(f,  z)  at  10,  22,  23.6,  25.1,  and  50  kilohertz.  The  10  kilohertz  component  attenu¬ 
ates  rapidly  with  distance,  whereas  components  at  the  lower  and  upper  transition 
frequencies  do  not  reach  steady  state  until  about  6  centimeters,  or  about  2.5  times 
the  mean  radius  of  the  tube.  These  figures  also  demonstrate  the  stable  spatial 
periodicity  of  the  high  frequency  spectral  components. 

3.9  Transformation  to  the  Time  Domain 

The  next  step  in  the  analysis  is  to  obtain  the  transient  response  of  the  tube 
as  a  function  of  distance  by  application  of  the  inverse  Fourier  transform  to  the 
transfer  functions,  w(f,z)  and  u(f}z).  Formally  this  is  defined  by: 

OO 

g(t,z)  =  2n  j  (3.50) 

—  OO 

We  will  perform  the  inversion  numerically  on  the  transfer  functions  by  means 
of  the  Inverse  Fast  Fourier  Transform  (IFFT),  which  is  a  widely  implemented 
version  of  the  discrete  Fourier  transform  (DFT).  The  theory  and  application  of 
numerical  versions  of  the  Fourier  transform  has  been  extensively  studied  by  many 
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authors,  including,  for  instance,  Bracewell  (1965),  Papoulis  (1977),  Burdic  (1984), 
Oppenheim  and  Schafer  (1989),  Jackson  (1986,  1991),  and  Strum  and  Kirk  (1995). 

Certain  restrictions  apply  when  using  numerical  tools  such  as  the  Fast  Fourier 
Transform  for  the  evaluation  of  integrals  such  as  equation  3.50.  First  of  all,  the 
DFT  uses  a  discrete  infinite  sum  as  an  approximation  to  the  continuous  function 
implied  in  3.50.  Thus  the  differential  df  in  the  integral  formulation  is  replaced 
by  an  increment  A /  in  the  discrete  version.  This  limits  the  resolution  which  may 
be  attained  with  the  discrete  approximation.  Nyquist  (1928)  proved  that  if  a 
continuous  time  function  g(t )  is  sampled  at  a  frequency  fs,  then  the  upper  limit  in 
frequency  space  for  spectral  components  of  the  Fourier  transform  of  g(t)  is  at  best 
equal  to  fa/2  ,  that  is,  one  half  of  the  sampling  frequency.  This  is  known  as  the 
Nyquist  frequency,  fN.  In  most  practical  applications,  the  bandwidth  is  restricted 
to  slightly  below  the  Nyquist  frequency.  The  infinite  DFT  may  be  viewed  as  a 
Fourier  series  representation  of  a  frequency  function  defined  over  a  finite  interval. 
This  series  will  also  represent  the  same  function  replicated  at  integral  multiples 
of  the  sampling  frequency.  If  the  original  time  function  is  not  band  limited  in 
frequency,  artificial  frequency  conponents  will  appear  in  the  DFT  spectrum.  This 
phenomenon  is  known  as  aliasing. 

A  second  constraint  results  as  a  consequence  of  replacing  the  infinite  sum  of 
the  DFT  with  a  finite  sum  of  N  terms  as  is  the  case  with  the  FFT.  The  spectrum 
obtained  with  the  finite  FFT  represents  a  function  which  is  replicated  in  the  time 
domain.  Thus  the  period  T  of  a  time  function  obtained  with  the  IFFT  is  limited 
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to  T  <  1/A/.  The  temporal  resolution  of  the  time  function  is  A t  =  T/N,  and 
therefore  Af  At  <  1/N. 

An  additional  constraint  applies  to  the  present  analysis  due  to  the  limita¬ 
tions  of  the  Donnell  shell  theory  as  described  in  the  analysis  of  the  thin  shell  in 
wavenumber-frequency  space,  that  is,  the  transfer  functions  w(f}z)  and  u(f,z ) 
must  be  band-limited  in  frequency  to  the  region  for  which  the  thin  shell  theory  is 
a  valid  model. 

We  recall  that  the  analysis  which  resulted  in  the  transfer  functions  for  w  and 
u  was  predicated  on  a  temporal  input  which  is  a  delta  function  of  time.  A  delta 
function  in  time  transforms  to  an  infinite  spectrum  in  frequency  space.  Before 
applying  the  inverse  transform  to  the  transfer  functions,  we  shall  band  limit  the 
spectrum  by  multiplying  each  transfer  function  by  an  appropriate  weighting  func¬ 
tion  H  (/)  chosen  so  as  to  assure  that  the  desired  constraints  are  realized.  Thus 
the  inversion  integral  3.50  will  take  the  form 

OO 

g(t,z)=  2a-  J  g(f,z)mf)e^’df  (3.51) 

— OO 

where  g(f,z )  represents  the  transfer  functions,  and  H(f)  is  a  weighting  func¬ 
tion  to  be  determined. 

3.9.1  Development  of  Weighting  Function 

The  characteristics  required  of  the  weighting  function  H(f)  are: 

(1)  A  uniform  ("flat”)  spectral  density  from  zero  to  some  upper  band  limit,  Fc 

(2)  A  rapid  and  smooth  transition  to  zero  spectral  density  for  frequencies  above 

Fc. 
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The  first  requirement  stems  from  the  need  to  simulate  the  spectrum  of  a  delta 
function  within  the  frequency  range  of  interest.  A  function  which  meets  this  re¬ 
quirement  is  the  rectangular,  or  ’’boxcar”  function,  which  is  defined  as: 


Hi(f)  = 


1  for  0  <  /  <  Fc 
0  for  /  >  Fc 


(3.52) 


Unfortunately  the  Fourier  transform  of  Hi(f)  does  not  meet  the  second  re¬ 
quirement.  The  abrupt  discontinuity  at  the  comer  frequency,  f  =  Fc,  will  produce 
the  so-called  ’’Gibbs’  Phenomenon”  (Strum  and  Kirk,  1995),  which  consists  of 
overshoot  and  ripples  in  the  transform  (time)  domain. 

One  function  which  does  meet  the  second  criterion  is  the  Gaussian  waveform 
(Burdic,  1984,  Bracewell,  1965).  The  Gaussian  function  may  be  expressed  as 

1  2 


I'hif)  =  exp  -a2 


(3.53) 


Equation  3.53  describes  a  Gaussian  frequency  distribution  centered  at  Fc  The 
parameter  a  controls  the  rate  at  which  this  function  falls  off  on  either  side  of  the 
cutoff  frequency  Fc. 

The  desired  weighting  function  is  synthesized  from  a  combination  of  3.52  and 
3.53: 

t 

Hi(f)  for  0  <  f  <  Fc 
H2(f)  for  f>Fc 
3.9.2  Weighted  Transfer  Functions 
Examples  of  the  weighted  transfer  functions  are  shown  in  figures  3-23  and  3-24 
for  axial  position  z  =  6  cm.  The  weighting  function  11(f  )  (multiplied  by  |  for 
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greater  visibility)  is  also  included  in  the  plot  of  w(f)  in  figure  3-23.  The  weighted 
transfer  function  is  labeled  if(/)  *  w(f)  in  this  figure.  We  note  that  the  weighted 
transfer  functions  make  a  gracefully  smooth  transition  to  zero  in  frequency  space. 
Thus  we  have  satisfied  the  requirement  discussed  earlier  that  the  spectrum  of  the 
transfer  functions  must  be  band-limited  before  inversion  to  the  time  domain  using 
the  IFFT.  Also  note  that  the  amplitude  of  u(f)  is  much  lower  than  the  amplitude 
of  w(f)  at  this  axial  position. 

3.9.3  Time  Domain  Response 

Application  of  the  IFFT  to  the  weighted  transfer  functions  shown  in  figures  3-23 
and  3—24  yields  the  radial  and  axial  displacements  as  functions  of  time  at  the 
selected  location,  z  =  6  cm.  This  result  is  shown  in  figures  3-25  and  3-26.  The 
effects  of  dispersion  on  the  propagation  of  waves  in  the  tube  is  clearly  evident 
in  these  plots.  The  high  frequency  components  of  the  transient  disturbance  arrive 
first.  These  are  followed  by  a  long  decaying  ’’tail”  in  which  the  frequency  gradually 
decreases.  The  frequency  components  in  this  ’’tail”  are  in  the  transition  zone 
discussed  earlier.  The  amplitude  of  the  radial  displacement  w(t)  is  about  50  times 
that  of  the  axial  displacement  u{t). 

The  resolution  in  the  time  domain  is  determined  by  the  number  of  points  N  and 
the  frequency  increment  A/,  as  mentioned  above.  The  frequency  increment  must 
be  large  enough  so  that  period  T  of  the  IFFT  includes  the  transient  of  interest. 
With  these  parameters  established,  the  minimum  number  of  points  Nm-m  depends 
on  the  desired  value  of  time  resolution,  At,  thus  Nm-m  =  1/(A /  At)  =  T/At. 
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3.10  Joint  Time-Frequency  Distributions 

The  fundamental  premise  which  supports  the  theory  of  Fourier  transforms  is 
that  the  process  which  generates  the  function  in  the  time  domain  is  stationary, 
that  is,  the  statistical  properties  of  the  process  are  invariant  in  time.  The  Fourier 
integral  itself  is  evaluated  formally  over  an  infinite  time  interval.  In  analyzing  real 
physical  data,  of  course,  the  sample  space  is  of  necessity  finite.  Real  signals  are 
limited  both  in  bandwidth  and  in  time.  However,  for  many  problems  of  interest, 
the  spectral  content  of  the  process  is  sufficiently  stable  over  the  time  span  of 
the  measurement  that  the  conditions  for  the  existence  and  validity  of  the  Fourier 
transform  are  satisfied. 

However,  there  are  many  processes  for  which  the  spectral  content  varies  sig¬ 
nificantly  during  the  time  sample  of  interest.  Common  examples  include  human 
speech,  bird  calls  and  other  biological  sounds,  music,  gravity  waves  in  large  bod¬ 
ies  of  water,  and  transients  in  dispersive  systems,  as  is  the  case  in  the  present 
study.  A  useful  construction  in  understanding  such  phenomena  is  the  concept  of 
a  joint  time-frequency  distribution,  in  which  the  energy  or  intensity  of  a  signal  is 
represented  simultaneously  in  time  and  in  frequency. 

The  topic  of  joint  time-frequency  distributions  has  received  intense  investiga¬ 
tion  by  the  signal  processing  community  in  recent  years.  An  excellent  review  of  the 
subject  has  been  published  by  Cohen  (1989).  A  study  of  joint  time-frequency  distri¬ 
butions  usually  begins  with  a  discussion  of  the  Wigner-Ville  Distribution  Function 
(WVDF).  Wigner  (1932)  devised  a  distribution  function  to  facilitate  his  study  of 
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quantum  statistical  mechanics.  Ville  (1948)  applied  the  same  type  of  distribution 
function  to  signal  analysis. 

Janse  and  Kaizer  (1983)  give  the  following  formulation  for  generalized  time- 
frequency  distributions: 

°° 

Cf{t,u,$)  =  ~  IJJ *u)$(S,T)f{u  +  ^)r(u- T-)dudrdi  (3.55) 

—  oo 

in  equation  3.55,  f(u)  is  the  time  signal,  f*(u)  is  the  complex  conjugate  of  the 
time  signal,  and  $(f,r)  is  a  kernel  function  which  determines  the  nature  of  the 
distribution  function.  The  WVDF  is  obtained  by  setting  $(£,t)  =  1  in  equa¬ 
tion  3.55.  Claasen  and  Mecklenbrauker  (1980)  established  that  all  time-frequency 
distribution  functions  may  be  obtained  by  convolution  of  the  Wigner  distribution 
with  some  window  function. 

In  the  present  study,  the  short  term  Fourier  transform,  or  spectogram,  will 
be  used  for  analysis  of  the  transient  response  of  the  tube  in  time-frequency  space. 
The  spectogram  is  a  well  known  tool  used  in  acoustical  studies  of  speech,  music, 
and  biological  events.  Authors  differ  in  their  definition  of  the  spectogram.  Cohen 
(1989),  Claasen  and  Mecklenbrauker  (1980),  and  Janse  and  Kaizer  (1983)  define 
the  spectogram  as  the  square  of  the  magnitude  of  the  short  term  Fourier  transform 
(SFT).  This  definition  emphasizes  the  power  spectral  density  of  the  signal.  Oppen- 
heim  and  Schafer  (1989)  and  Rabiner  and  Schafer  (1978)  define  the  spectogram  as 
the  magnitude  of  the  SFT.  We  shall  use  the  latter  definition  in  the  present  study 

The  short-time  Fourier  transform  is  given  by  (Claasen  and  Mecklenbraucker, 
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1980): 


/oo 

e'"T  f(r)  h(r  -  t)  dr 

•OO 


(3.56) 


in  which  h{r  —  t)  is  a  window  function  centered  at  time  t,  and  /(r)  is  the 
signal.  Then  the  spectogram  is  Sf(t,u)  =  \ Ft  (ur) j .  In  this  study  we  shall  employ 
an  algorithm  for  the  spectogram  which  is  implemented  in  the  MATLAB  (1994) 
system  for  numerical  analysis.  The  spectogram  is  developed  by  sliding  a  window 
along  the  time  axis.  At  each  discrete  position  of  the  window,  the  SFT  is  com¬ 
puted.  The  result  is  a  matrix  in  which  each  column  represents  the  spectrum  at  a 
particular  instant  of  time.  The  rows  in  the  matrix  contain  the  complex  values  of 
the  spectogram  at  specific  frequencies. 


The  Uncertainty  Problem  There  is  an  inherent  fundamental  limitation  in  the  ap¬ 
plication  of  time-frequency  distributions,  namely,  the  duration  of  a  signal  and 
its  bandwidth  are  reciprocally  related.  This  is  often  referred  to  the  ’’uncertainty 
principle,”  in  analogy  to  the  famous  Heisenberg  uncertainty  priciple  in  quantum 
mechanics.  As  generally  expressed,  the  uncertainty  principle  states  that  the  dura¬ 
tion  of  a  signal  and  the  spectral  bandwidth  of  its  Fourier  transform  are  bounded 

by 

(AT)  (AC)  >  -t  (3.57) 

where  AT  is  the  length  of  the  signal  in  seconds,  and  AF  is  the  spectral  band¬ 
width  in  hertz. 

Proof  of  this  may  be  found  in  Papoulis  (1977),  Claasen  and  Mecklenbrauker 
(1984)  and  Elmore  and  Heald  (1985).  Cohen  (1989)  offers  the  following  quotation 


from  Skolnik  (1980):  ”The  use  of  the  word  ’uncertainty’  is  a  misnomer,  for  there 
is  nothing  uncertain  about  the  ’uncertainty  relation.’.  .  .  It  states  the  well-known 
mathematical  fact  that  a  narrow  waveform  yields  a  wide  spectrum  and  a  wide 
waveform  yields  a  narrow  spectrum  and  both  the  time  waveform  and  the  frequency 
spectrum  cannot  be  made  arbitrarily  small  simultaneously.”  Applications  of  joint 
time-frequency  distributions  are  described  by  Yen  (1987),  Janse  and  Kaizer  (1983), 
Rossano,  Shin,  and  Hamilton  (1990),  and  Rao,  Taylor,  and  Harrison  (1990). 

The  uncertainty  relation  enters  into  the  computation  of  the  time  domain  re¬ 
sponse  in  terms  of  the  cutoff  frequency  of  the  weighting  function,  H(f).  If  the  ef¬ 
fective  bandwidth  of  the  weighting  function  is  Fc,  then  the  minimum  signal  length 
in  the  time  domain  is  given  by  AT  =  1/(47T  Fc). 

Spectogram  Parameters:  We  now  consider  some  of  the  factors  which  enter  into  the 
computation  of  the  spectogram  of  discrete,  finite  time  series.  In  the  experimental 
study  of  time-dependent  phenomena,  measurements  may  be  made  as  a  continuous 
function  of  time  with  an  analog  device,  such  as  a  conventional  tape  recorder,  or 
as  a  discrete  function  of  time  using  a  digital  device  such  as  an  analog  to  digital 
(A/D)  converter.  In  the  first  case,  the  temporal  resolution  is  determined  by  the 
bandwidth  of  the  measurement  device.  In  the  digital  case,  the  limitations  are  set 
by  the  sampling  frequency  Fs  (Hertz)  and  the  maximum  number  of  data  points 
which  may  be  recorded  with  the  instrument,  N.  The  bandwidth  in  the  discrete  case 
must  be  limited  by  anti-aliasing  filters  to  a  value  less  than  the  Nyquist  frequency, 
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The  sampling  interval  At  is  the  reciprocal  of  the  sampling  frequency,  At  = 
1/FS,  and  the  record  length,  or  maximum  signal  duration,  is  Tmax  =  N  •  At.  (But 
note  that  for  a  signal  that  is  bounded  in  time,  such  as  a  transient,  the  effective 
signal  duration  may  be  artificially  extended  by  adding  zeros  to  the  time  series, 
a  technique  known  as  ’’zero  padding”).  The  uncertainty  relation  states  that  the 
minimum  bandwidth  of  the  Fourier  transform  is 


AF  > - 

47riV  At 


Fs 

4itN 


There  is  an  inherent  trade-off  between  the  frequency  and  time  resolution  in  a 
spectogram.  The  frequency  resolution  Af  of  the  spectrum  of  the  entire  time  series 
of  N  points  is  given  by  A /  =  =  jf-  However,  in  the  algorithm  employed  for 

the  computation  of  the  spectograms,  the  original  time  series  of  N  data  points  is 
divided  into  sections  of  length  N/ft,  and  then  each  section  is  multiplied  by  a  window 
function.  (The  Hanning  window  is  used  in  this  study.)  The  frequency  resolution  of 
the  spectogram  is  determined  by  the  length  of  the  original  sampling  interval  and 
the  length  of  the  section,  thus  A fsp  =  1  /Njft  At.  Thus  the  frequency  resolution 
is  reduced  in  the  spectogram  in  the  ratio  N/Nfft.  Put  another  way,  spectral  fines 
which  were  narrow  in  the  spectrum  of  the  whole  time  series  will  be  broadened 
or  smeared  in  the  spectogram.  On  the  other  hand,  the  time  domain  resolution  is 
approximately  determined  by  the  same  ratio,  N/Nfft  which  determines  the  number 
of  points  in  time  at  which  the  FFT  is  computed.  This  problem  is  endemic  to  all 
joint  time-frequency  distributions  (Classen  and  Mecklenbrauker,  1980, 1984).  Thus 
a  judicious  compromise  must  be  made  in  selecting  Nfft  for  the  spectogram. 
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Another  factor  which  influences  the  appearance  of  spectograms  is  known  as 
’’overlap.”  It  is  customary  to  include  a  certain  amount  of  overlap  between  con¬ 
tiguous  sections  of  the  subdivided  time  series.  If  the  number  of  overlap  samples 
is  N0i ,  then  the  number  of  discrete  time  points  in  the  spectogram  is  given  by  the 
integer  part  of  (N  —  Nol)/(Nfft  —  Noi).  Generally  the  overlap  is  set  to  Nfft/2,  so 
that  half  of  the  samples  of  the  J  th  section  are  included  in  the  ( J  +1)  th  section. 
For  this  case  the  number  of  time  points  is  the  integer  part  of  —  1. 

Spectograms  of  the  waveforms  displayed  in  figures  3-25  and  3-26  are  presented 
as  contour  plots  in  figures  3-27  and  3-28  respectively.  In  this  representation,  the 
absolute  value  of  the  FFT  results  in  the  spectogram  matrix  are  normalized  to 
the  global  maximum.  Contour  lines  connect  points  of  equal  magnitude.  The 
spectograms  are  shaped  somewhat  like  a  boomerang,  with  a  steep  ridge  parallel  to 
the  frequency  axis  occurrring  early  in  time.  The  long  tail  parallel  to  the  time  axis 
represents  the  low  frequency  components.  The  spectogram  results  agree  with  the 
description  of  these  waveforms  earlier  in  the  section  ’’Time  Domain  Response.” 
The  longer  tail  in  the  spectogram  of  the  axial  component  reflects  the  increased 
ralative  amplitude  of  the  low  frequencies  evident  in  figure  3-26. 

Figures  3-29  and  3-30  are  three  dimensional  contour  plots  of  the  spectograms. 
In  the  three  dimensional  plots,  contour  planes  of  equal  amplitude  are  drawn  parallel 
to  the  time-frequency  plane.  The  dispersive  nature  of  the  wave  propagation  is 


clearly  evident  in  the  contour  plots. 
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CHAPTER  4 
Composite  Systems 

4. 1  Introduction 

The  models  developed  in  Chapters  2  and  3  are  brought  together  here  to  for¬ 
mulate  the  description  of  various  composite  systems  of  interest.  This  will  include 
models  comprised  of  all  fluid  components,  and  systems  in  which  elastic  media  are 
incorporated.  The  transfer  matrix  approach  described  in  Chapter  1  will  be  used 
for  this  study. 

4.2  Fluid  Systems 

4.2.1  Parameters  for  Two  Fluid  Systems:The  behavior  of  a  system  which  includes 
two  different  fluids  will  depend  on  the  density  and  sound  speed  in  the  fluids  as  well 
as  the  frequency  and  wavenumber  of  the  excitation.  Analysis  of  such  systems  is 
simplified  by  defining  several  non-dimensional  ratios  of  various  physical  param¬ 
eters.  For  a  two  fluid  system,  the  subscript  1  will  denote  the  inner  fluid,  and 
subscript  2  will  represent  the  outer  fluid.  We  now  define  the  following  ratios: 

Ratio  of  sound  speeds: 

Rc  =  — 

Cl 


Ratio  of  densities: 
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Ratio  of  specific  acoustic  impedances: 


Rz 


P2.^2  _  D  I? 

-  lhC-Llp 

Pi  Cl 


Table  4-1,  Parameter  Matrix  for  Two-Fluid  Problems 


Sound  Speed 

Ci  <  c2 

Cl  =  c2 

Cl  >  c2 

Wavenumber 

ki  > 

ki  =  k2 

Cl  >  c2 

Density 

Pi  <  P2 

RZ,RP,  Rc  >  1 

RZ,RP  >  1,.RC  =  1 

Rz  =?,RP  >1,RC<1 

Pi  =  P2 

Rz  >  1,  Rp  =  1,RC  >  1 

Rz,Rp,Rc  =  1 

Rz  <  1,  Rp  =  1,  Rc  <  1 

Pi  >  P2 

Rz  =?,RP  <  1,-Rc  >  1 

Rz,Rp<l,Rc  =  l 

Rz  >  1  ,RP,RC  <  1 

Table  4-1,  ’’Parameter  Matrix  for  Two-Fluid  Problems”,  indicates  that  there 
are  nine  combinations  of  fluid  density  and  sound  speed  which  must  be  considered 
in  two  fluid  problems.  The  relative  magnitudes  of  the  three  ratios  defined  above 
are  also  indicated  in  Table  4-1. 

There  will  be  two  acoustic  wavenumbers,  ki  =  u)/c\,  and  k2  =  w/c2,  and  two 
radial  wavenumbers,  A&  =  k\  -  k2,  and  k2r2  =  k\-  k2.  When  the  frequency  and 
axial  wavenumber  are  both  independent,  we  may  define  a  variable  cz  =  lo / kz  which 
represents  the  speed  at  which  a  wave  propagates  in  the  z  direction.  (This  is  some¬ 
times  called  the  ’’trace  velocity”).  We  now  define  the  ratio  of  radial  wavenumbers 
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as: 


For  propagation  to  occur  in  both  fluids,  kT 2  and  kr i  must  be  real.  This  will  be 
true  if  both  C\  and  c2  are  less  than  cz.  For  a  system  driven  at  the  inner  boundary 
of  fluid  1,  if  ci  <  cz  and  c 2  >  cz,  propagation  will  occur  in  fluid  1,  but  not  in  fluid 
2.  If  Ci  >  cz,  propagation  will  not  occur  in  either  fluid. 

4.2.2  Fluid  Tube  in  External  Fluid:The  first  case  consists  of  a  fluid  tube  driven 
at  its  inner  surface  and  imbedded  in  an  external  fluid.  The  external  medium  has 
a  real  component  of  radiation  impedance,  enabling  the  radiation  of  energy. 

Consider,  then,  a  fluid  tube  of  inner  radius  and  outer  radius  a2,  with  fluid 
density  pi,  sound  speed  Ci,  and  wavenumber  k\  =  — .  The  fluid  tube  is  in  intimate 
contact  with  the  external  fluid  of  density  p2  and  sound  speed  c2  at  the  boundary, 
r  =  a2.  The  wavenumber  in  the  external  fluid  is  fc2  =  The  input  to  the  system 
is  a  traveling  wave  of  pressure  or  velocity  of  the  form  PQ  exp i(kzz  —  cot),  where 
P0  is  an  amplitude  coefficient,  kz  is  an  arbitrary  axial  wavenumber,  and  u ;  is  the 
angular  frequency. 

This  model  could  represent,  for  example,  a  cylindrical  transducer  or  towed 
array  in  which  a  ceramic  cylinder  or  array  of  cylinders  is  enclosed  in  a  soft  elastomer 
(with  a  negligible  shear  modulus),  immersed  in  sea  water. 

The  boundary  conditions  which  must  be  satisfied  at  the  interface  r  =  a2  are: 
(1)  continuity  of  the  state  variables  of  pressure  and  radial  velocity,  and  (2)  waves 


in  the  fluid  tube  and  in  the  external  fluid  must  have  the  same  axial  wavenumber 
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kz  if  propagation  is  to  occur. 

The  pressure  and  radial  velocity  at  r  —  a2  are  related  by  the  spectral  impedance 
of  the  external  fluid;  that  is,  Pa2  =  Z2Va2,  where 

^  _  -ip2C2fc2go1)(fcr2Q2) 

kT2H^\kT2a2) 

Using  the  equations  for  a  fluid  tube, 


Pa* 

K2 


=  [Mab] 


Pa , 


where 


[Mab]  = 


Ay  By 
Cl  Dy 


and  Ay.By,CyDy  are  the  ABCD  parameters  for  the  fluid  tube  as  defined  in 


Chapter  2. 


Solving  for  Vai  we  find 


Vai  =  Pc 


(Ay  -  CyZj) 
(DyZ2  -  By) 


and  the  pressure  at  r  =  a2  is 


Pa 2  —  Pc 


oi 


Ay  + 


By  (Ay  —  C1Z2) 
(DyZ2~By) 


thus,  all  of  the  state  variables  are  determined  from  either  the  applied  pressure 
Pai  or  velocity  Vai ,  the  fluid  tube  dimensions  ay  and  a2,  and  the  acoustic  properties 
of  both  media. 

The  ABCD  parameters  for  the  fluid  tube  are  given  by 


7 r. 


Ay  —  ^  kri  ayQyo  (kTi ) 
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B\  —  i  ^  for i  aiQoo  (&ri  )  ^1 

c,=if%fQn(iv,) 

7T 

Z?1  =  — fcri &1Q10  (^ri  ) 


where  Zi  = 


picik} 

kr^ 


Now  consider  the  spectral  impedance  at  boundary  d\, 


7  _  P*i 

Za'~v7 


ai 


D\Z%  —  Bi 

A  -  C^z2 


Substituting  the  explicit  expressions  into  4.1  and  simplifying  yields 


Za j  —  %Z\ 


r2 


For  convenience  we  define  the  following  ratios: 


Ratio  of  Hankel  functions: 


(4.1) 


(4.2) 


Hrat(x)  = 


flf(*) 


Ratio  of  radii: 


ai 

Note  that  =  Rc,  and  also  2^L  = 


Substituting  these  quantities  into  4.2  yields 


Zai 


=  iZi 


[jRp^Qoi  (KJ  Hrat  ( kT2sa )  —  Qoo  (fcn) 
^Qio(^ri)  (^ri )  Hrat  (kr2sa) 


(4.3) 


The  response  of  this  system  to  a  specific  excitation  will  depend  on  the  wavenum¬ 
ber  and  frequency  of  the  excitation  [kz  and  u>),  as  well  as  the  acoustical  properties 
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of  the  fluid  tube  and  the  external  medium.  Several  limiting  cases  of  interest  will 
be  studied  using  equation  4.3. 


Case  I  -  Zero  Axial  Wavenumber  The  simplest  case  is  that  for  which  there  is  no 
axial  wave  propagation,  that  is,  kz  =  0.  The  excitation  at  r  =  d\  is  a  pressure  or 
velocity  which  is  uniform  in  amplitude  along  the  z  axis,  varying  harmonically  at 
frequency  u. 

For  this  case,  kri  =  k\,  kr2  =  k%,  and  =  Rc-  also,  Z\  =  p\C\.  The  spectral 

impedance  at  r  =  a*  (normalized  to  the  specific  acoustic  impedance  of  medium  1 , 
PiCj),  is  found  from  4.3  as 


Ja  i 


RzQoi  (h)  Hrat  (^)-Qoo(h) 


(4.4) 


Pid  [Qio  (fci)  -  RzQu  (ki)  Hrat  (^)] 

The  variables  for  this  special  case  are  the  non-dimensional  wavenumber  in  the 
fluid  tube,  Aqai,  the  ratio  of  specific  acoustic  impedances  Rz,  the  ratio  of  sound 
speeds  Rc,  and  the  ratio  of  radii,  s. 

Now  for  large  values  of  the  arguments  in  equation  4.4,  that  is,  k\a\,  k^a^,  and 
fc2a.2  are  all  greater  than  2tt,  substitution  of  the  asymptotic  expressions  developed 
in  Chapter  2  leads  to  the  following  approximation: 


ljm  (  Za^  \  —  itan[fci(a2  ~  Qi)] 

fcioi-^oo  l  p\C\)  1  —  iRz  tan[fci(a2  —  Qi)] 


(4.5) 


With  these  simphfications,  the  spectral  impedance  at  the  inner  surface  of  the 
fluid  tube  is  determined  by  Rz ,  the  ratio  of  acoustic  impedances,  and  the  argument 
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of  the  tangent  function,  fci(a2  —  ax),  the  product  of  the  wavenumber  in  fluid  1  and 
the  thickness  of  the  tube.  Two  limiting  cases  of  interest  will  be  discussed. 

(1)  Tube  thickness  an  integral  number  of  half- wavelengths:  If  the  tube  thickness 
02  —  ax  is  an  integral  number  of  half-wavelengths,  that  is, 

nit  „  „  „ 

a2  -  Oi  =  — ,  n  =  1,2,3. . . 

kj 

then  Zai  — >  p2c2;  the  characteristic  impedance  of  the  external  fluid. 


(2)  Tube  thickness  an  odd  number  of  quarter-wavelengths:  If 


q2  —  Oj 


(2n  -1-  l)7r 
4fci  ’ 


n  =  1,2,3. . . 


(4.6) 


then  Zai  — *  •  This  relationship  suggests  an  interesting  method 


for  ’’matching"  the  impedance  of  a  source,  such  as  a  ceramic  ring,  to  the  specific 
acoustic  impedance  of  the  external  fluid.  If  Zai  represents  the  impedance  of  the 
source,  the  fluid  in  the  tube  would  be  chosen  such  that 


PlCl 


\J  Zai 


P2C2 


Of  course,  the  assumptions  made  in  obtaining  equation  4.5  must  hold,  and 
exact  ’’matching"  will  occur  only  at  frequencies  satisfying  equation  4.6.  For  other 
combinations  of  frequency  and  tube  thickness,  the  impedance  Zai  will  be  complex. 
In  rectangular  components,  this  may  be  expressed  as 

Vm  {  Zaj_\  _  Rz  +  tan2  [Aq  (q2  -  tan  [fcx  (a2  -  aQ] 

°i  — >00  ypjCi  J  1  +  R%  tan2  [kx  (a2  —  ai)]  l  +  R%  tan2  \kx  (a2  —  ax)] 
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Returning  to  the  general  case  described  by  equation  4.4,  we  will  discuss  the 
behavior  of  the  normalized  specific  acoustic  impedance  as  a  function  of  the  non- 
dimensional  wavenumber  k\a,\  for  various  combinations  of  the  parameters  RC,RZ, 
and  s. 

(1)  Variation  with  tube  thickness,  ax(s  -  1):  For  this  case,  we  assume  a  sound 
speed  ratio  of  Rc  =  2,  and  an  impedance  ratio  Rz  =  2.  Figure  4-1  presents  the 
variation  in  impedance  as  a  function  of  k\a\  for  three  values  of  the  radius  ratio, 
s  =  1,  2,  and  4.  The  magnitude,  phase  angle,  and  rectangular  components  of  the 
impedance  are  shown. 

The  degenerate  case  5  =  1,  shown  as  a  solid  line,  represents  the  impedance 
of  the  external  medium  alone,  as  discussed  in  Chapter  2.  (For  this  case,  there 
is  no  fluid  tube  in  the  system).  For  large  values  of  the  argument,  the  spectral 
impedance  approaches  the  specific  acoustic  impedance  of  the  external  medium, 
p2c2.  For  s  =  2  (dotted  line)  and  s  =  4  (dashed  line)  the  impedance  functions 
exhibit  cyclic  behavior  as  a  function  of  acoustic  wavenumber  (or  frequency).  The 
period  of  the  fluctuations  decreases  with  increasing  tube  thickness,  in  a  fashion 
similar  to  that  discussed  for  the  fluid  tube  in  Chapter  2.  The  impedance  magnitude 
varies  from  Rz  to  l/Rz,  as  in  the  simplified  case  described  above.  The  imaginary 
component  of  impedance  is  skewed  in  a  triangular  fashion  The  real  component 
waveform  is  less  than  unity  for  a  greater  part  of  the  cycle,  and  exceeds  unity 
for  the  remainder  of  the  cycle.  At  large  values  of  the  argument,  a  steady-state 
condition  is  reached  for  the  amplitude  of  the  impedance  variations. 
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(2)  Variation  with  acoustic  impedance  ratio,  Rz :  In  figure  4-2,  a  sound  speed 
variation  of  two  to  one  is  assumed  ( Rc  =  2),  and  the  fluid  tube  thickness  is  twice 
the  inner  radius  (s  =  2).  The  spectral  impedance  is  shown  for  three  values  of  the 
impedance  ratio,  Rz  =  |,1,  and  2.  Since  the  sound  speed  ratio  is  fixed  at  2,  the 
impedance  ratio  Rz  is  determined  by  the  fluid  density  in  each  medium.  Thus  for 
Rz  =  1,  p2  —  etc. 

Again,  a  cyclic  variation  of  spectral  impedance  with  wavenumber  is  observed. 
The  phase  of  the  variation  is  determined  by  the  magnitude  of  the  impedance 
ratio.  (Compare,  for  example,  the  solid  curves  for  Rz  =  \  with  the  dashed  curves 
for  Rz  =  2.  )  The  magnitude  of  the  variation  also  depends  on  the  value  of  Rz, 
being  lowest  for  Rz  —  1,  and  increasing  for  both  larger  and  smaller  values  of 
the  impedance  ratio.  The  same  comments  regarding  the  waveforms  of  the  real  and 
imaginary  components  as  made  in  the  preceding  case  apply  to  the  present  analysis. 

(3)  Variation  with  sound  speed  ratio,  Rc:  For  this  example,  the  impedance  ratio 

is  fixed  at  1.5,  and  the  ratio  of  radii  is  again  2.  Figure  4-3  presents  the  spectral 
impedance  versus  frequency  for  Rc  =  1,  and  10.  When  the  sound  speed  in  the 

external  fluid  exceeds  that  in  the  fluid  tube  (dashed  line),  the  impedance  variation 
is  greatest  at  small  values  of  Aqaj,  decreasing  with  increasing  wavenumber  until 
a  steady  state  is  reached.  The  reverse  is  true  for  Ci  >  c%,  where  the  impedance 
variation  builds  up  to  a  steady  state  value  as  the  frequency  increases.  The  same 
asymmetry  previously  noted  in  the  real  and  imaginary  components  is  apparent  for 
this  case. 
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In  summary,  there  are  three  parameters  which  affect  the  normalized  spectral 
impedance  of  the  fluid  tube  -  external  fluid  system.  The  tube  radius  ratio  5  controls 
the  periodicity  of  the  cyclic  impedance  variations,  with  the  period  decreasing  as 
s  increases.  The  fluid  impedance  ratio  Rz  and  sound  speed  ratio  Rc  primarily 
influence  the  magnitude  of  the  cyclic  impedance  variations.  In  all  cases,  the  cyclic 
variations  reach  a  steady  state  for  k\a\  27T. 

Case  II  -  Finite  Axial  Wavenumber:  The  more  general  case  in  which  the  axial 
wavenumber  is  finite  will  include  both  propagating  waves  as  described  earlier,  and 
standing  waves  in  which  the  excitation  may  be  of  the  form  P0(cos  kzz)(cosujt).  The 
conditions  for  the  existence  of  a  standing  or  a  propagating  wave  in  the  fluid  tube 
were  described  earlier,  that  is,  k\  >  kz.  The  radial  spectral  impedance  at  r  =  d\ 
will  be  given  by  equation  4.3.  The  behavior  of  the  impedance  will  be  similar  to  that 
described  for  Case  I  when  both  radial  wavenumbers  kr\  and  kr%  are  real.  When 
ki  <  kz  an  evanescent  wave  will  occupy  the  fluid  tube.  This  may  ’’leak”  into  the 
external  medium  if  the  fluid  tube  is  sufficiently  thin. 

4.2.3  Multiple  Concentric  Fluid  Tubes: We  now  consider  a  system  comprised  of 
any  number  N  of  concentric  fluid  tubes  immersed  in  an  external  fluid,  as  illustrated 
in  Figure  1-1.  The  first  fluid  tube  extends  radially  from  r  =  a\  to  r  =  a2  and 
contains  a  fluid  of  density  p\  and  sound  speed  C\.  The  j  th  tube  is  bounded  by  the 
inner  radius  CLj  and  outer  radius  a^+i,  with  fluid  of  characteristic  impedance  pjCj. 
The  outermost  fluid  tube  extends  from  radius  a^-i  to  a at,  where  a at  is  the  radius 
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of  the  interface  between  the  N  th  tube  and  the  external  fluid. 

Now  let  Pi  and  Vi  be  the  state  variables  of  pressure  and  radial  velocity  applied 
to  the  innermost  fluid  tube  at  r  =  a\.  The  pressure  and  radial  velocity  at  the 
interface  r  =  a2  are  given  by 


P‘1 

A! 

Bi 

Pi 

V*  _ 

Q 

Di 

.  . 

(4.7) 

The  values  of  P2  and  V2  from  equation  4.7  become  the  input  for  the  second 
fluid  tube  of  inner  radius  a^.  so  that  the  transfer  function  from  ai  to  <2.3  may  be 
expressed  as 


P3 

A2 

b2 

A\ 

Bi 

Pi 

V3 

C2 

d2 

Cl 

Di  _ 

Vi 

This  process  is  repeated  for  each  fluid  tube  in  the  system,  with  the  result 


Pn 

Ax 

Bx 

Pi 

VN 

Cx 

Dx 

Vi 

where  P n  and  Vn  are  the  state  variables  at  r  =  a^,  the  interface  with  the 
external  fluid,  and  Pjv  =  VnZn,  in  which  Zn  is  the  spectral  impedance  of  the 
external  fluid,  equation  2.57. 

The  overall  transfer  matrix  is  obtained  by  cascading  the  matrices  for  each 
component,  viz 


Ax  Bx 

An-  1  Bn-  1 

r  -| 

_ 1 

r  -1 

Ax  Bx 

1 

Cn-i  Dn~  1 

- 

Cj  D i 

- 

- , 

Q 

o' 

_ 1 

Care  must  be  taken  when  using  this  model  that  the  correct  ABCD  expressions 
are  used  in  each  component,  depending  on  whether  the  radial  wavenumber  k3  is 
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real  or  imaginary.  The  conditions  discussed  earlier  in  the  section  ’’Parameters 
for  Two  Fluid  Systems”  must  be  evaluated  at  each  radial  interface  to  determine 
whether  propagating  or  evanescent  waves  will  exist  in  the  succeeding  fluid  tube. 

4.2.4  Systems  Including  a  Fluid  Cylinder:When  a  composite  cylindrical  system 
includes  a  fluid  cylinder,  this  component  may  exert  a  major  influence  on  wave 
propagation.  When  the  radial  wavenumber  in  the  fluid  cylinder  is  real,  equation 
2.14  and  Figure  2.3  show  that  the  spectral  impedance  of  the  cylinder  will  cycle 
from  zero  to  infinity  as  a  function  of  the  argument  krca  where  krc  is  the  radial 
wavenumber  and  a  is  the  radius  of  the  cylinder.  This  establishes  the  range  of 
boundary  conditions  at  the  cylinder  surface,  ranging  from  zero  velocity  (’’blocked” 
pressure)  to  zero  pressure  (’’free”  surface). 

When  the  axial  wavenumber  kz  is  much  greater  than  the  acoustic  wavenumber 
in  the  fluid  cylinder,  equation  2.17  predicts  a  mass-like  spectral  impedance  which 
increases  with  frequency  and  decreases  with  axial  wavenumber. 

Fluid  Cylinder  within  Fluid  Tube:  We  now  consider  a  fluid  cylinder  of  radius  ai 
and  specific  acoustic  impedance  p\C\  which  is  concentric  with  a  fluid  tube  of  outer 
radius  a2  and  filled  with  a  fluid  of  impedance  p2c2.  This  system  is  driven  by  a 
wave  of  frequency  u  and  wavenumber  kz  applied  at  the  outer  radius.  The  transfer 
matrix  for  the  fluid  tube  is 


Pl 

a.21 

821 

P2 

_  Vi  _ 

_C21 

821 

v2  _ 

where  the  ABCD  functions  are  as  defined  by  equations  2.37  to  2.40.  The 
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pressure  and  velocity  at  r  =  ai  are  related  by  P\jV\  —  Zjc,  the  spectral  impedance 
of  the  fluid  cylinder.  Solving  for  the  spectral  impedance  at  the  input,  r  =  a2,  we 


have 

£2  _  2  —  ~  ^21 

^2  2  -^21  _  CziZfc 

When  the  impedance  of  the  fluid  cylinder  is  zero,  Z2  is  given  by 


(4.8) 


„  |  B21  .P2UJ  [Jo  {kr 2O2)  ho  (kr2al)  —  Jo  {kr2al)  ^0  (kr2a2)] 

//2  ^  ¥  Q  - zzz  — ^ - - - - - — .  — . - 

A21  kr 2  [Jx  (kr2a2)  ho  (K2CL1)  —  Jo  (kr 2ai)  hi  (fc^a^)] 
Similarly,  when  the  impedance  of  the  fluid  cylinder  is  infinite,  Z2  is 


£21  _  iP2^_  [Ji  (fcr2fli)  hp  (fcr2a2)  -  Jo  (£202)  hi  (fer2ai)] 

C2I  fcr2  [Jl  (fcr202)  hi  (fcr2ai)  -  Ji  (fc^Ol)  hi  (fc^)] 


The  impedance  Z2  is  always  reactive,  since  there  is  no  dissipation  or  energy 


radiation  in  the  system.  The  magnitude  and  phase  of  the  spectral  impedance  will 


vary  with  the  dimensions  ai  and  02,  as  well  as  the  frequency  and  wavenumber  of 
the  excitation. 


4.3  Fluid  Loaded  Elastic  Shells 

One  of  the  most  widely  studied  topics  in  the  area  of  structural  acoustics  is  that 
of  fluid  loaded  elastic  shells.  The  preponderance  of  work  has  been  concerned  with 
cylindrical  elastic  shells,  first,  because  of  the  multitude  of  practical  problems  which 
may  be  described  by  this  configuration  (several  examples  were  cited  in  Chapter  1), 
and  second,  because  the  separability  of  the  wave  equation  in  cylindrical  coordinates 
renders  such  problems  analytically  tractable. 

Junger  and  Feit  (1986)  devote  substantial  attention  to  both  finite  and  infinite 
fluid  loaded  cylindrical  shells,  including  sound  radiation  by  driven  shells  and  the 
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scattering  of  incident  sound  by  elastic  cylinders.  Their  book  ’’Sound,  Structures 
and  Their  Interaction”  is  an  excellent  starting  point  for  the  study  of  this  topic. 
Most  authors  use  thin  shell  theories,  such  as  described  in  Chapter  3,  to  repre¬ 
sent  the  elastic  cylinder  in  their  analyses  of  fluid  loaded  shells.  One  exception  is 
Greenspan  (1960)  who  used  elasticity  theory  in  his  study  of  an  elastic  shell  with 
fluid  loading  on  both  surfaces  and  driven  by  a  spatially  arbitrary,  time  harmonic 
pressure  applied  to  the  inner  surface.  A  more  recent  example  is  Pathak  (1993)  who 
applied  elasticity  theory  to  study  the  radiated  pressure  from  an  infinite  length  thick 
cylindrical  shell  harmonically  driven  by  a  ring  force.  The  present  discussion  will 
employ  the  axisymmetric  thin  shell  model  developed  in  Chapter  3. 

4.3.1  Elastic  Shell  in  Fluid  System:Now  consider  a  consider  a  system  with  a  fluid 
tube  of  inner  radius  a*  and  outer  radius  a 2,  with  an  elastic  shell  at  this  outer 
radius.  Since  the  thickness  of  the  elastic  shell  is  very  small  compared  to  its  radius, 
we  may  neglect  it  when  computing  fluid  impedances.  A  second  fluid  tube  extends 
radially  from  a 2  to  a 3  where  it  meets  the  external  fluid.  The  specific  acoustic 
impedance  of  the  first  fluid  is  P1C1,  the  second  fluid  is  P2C25  and  the  external  fluid 
is  P3C3.  This  model  could  represent,  for  example,  a  sonar  transducer  comprising  a 
piezoelectric  ceramic  cylinder  surrounded  by  a  coupling  fluid,  such  as  castor  oil  or 
isoper,  a  thin  metal  casing,  and  a  protective  outer  elastomeric  jacket  of  a  material 
with  negligible  shear  modulus.  The  external  fluid  in  this  case  would  be  water.  The 
system  input  at  r  =  a\  would  be  provided  by  the  piezoelectric  ceramic  cylinder. 


Following  the  method  described  under  ’’Multiple  Concentric  Fluid  Tubes,”  the 
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model  for  such  a  system  is: 


Pa 
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Bz 

1 

Ax 

Pi 

Pi 

V, 

_c3 

Dz 

0  1 

Cx 
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in  which  Pi  and  V1  are  the  pressure  and  radial  velocity  at  the  interface  between 
the  ceramic  cylinder  and  the  first  fluid  tube;  A\  B\  C\  D\  represent  the  first  fluid 
tube,  A3  Bo  Cz  Dz  represent  the  second  fluid  tube,  and  Zs  is  the  spectral  impedance 
of  the  shell,  from  Eq.  3.8.  The  pressure  and  radial  velocity  at  r  =  a3  are  related 
by  the  spectral  impedance  of  the  external  fluid,  P4  =  Z4V4,  where  Z\  is  computed 
by  substitution  of  the  characteristics  of  the  external  fluid  into  Eq.  2.57.  Thus  the 
state  and  field  variables  of  the  entire  system  may  be  determined. 

A  generalized  discussion  of  this  system  is  not  feasible  due  to  the  large  number 
of  parameters  involved  (three  radii,  three  fluids,  and  the  shell  characteristics). 
However,  a  specific  design  problem  is  easily  approached  following  the  methods 
described  in  this  and  earlier  chapters. 

4.3.2  Elastic  Shell  With  Fluid  Cylinder:The  problem  of  wave  propagation  in  a 
fluid  filled  elastic  shell  has  been  approached  by  many  authors,  including  Thomson 
(1953),  Skalak  (1956),  White  and  Sawley  (1972),  Kumar  (1966,1972).  Merkulov, 
Prikhod’ko  and  Tyutekin  (1978),  Fuller  (1981,  1988),  Ujihashi  et  al  (1987),  Fuller 
and  Fahy  (1982),  and  more  recently,  Brevart  and  Fuller  (1994).  The  great  interest 
in  this  problem  stems  from  practical  considerations  such  as  noise  (e.g.  "water 
hammer” )  in  piping  systems  and  the  reduction  of  noise  from  the  exhaust  of  internal 


combustion  engines. 


Ill 


The  propagation  of  free  waves  in  a  fluid  filled  thin  elastic  cylinder  was  studied 
by  Fuller  and  Fahy.  They  computed  the  dispersion  curves  for  circumferential  modes 
of  order  n  =  0,1  and  2  for  steel  and  hard  rubber  shells  filled  with  water.  The 
coupling  of  the  fluid  to  the  shell  complicates  the  solution  for  the  complex  roots  of 
the  characteristic  equation  because  the  roots  are  embedded  in  the  argument  of  the 
Bessel  functions  which  represent  the  fluid  loading  term  (that  is,  the  argument  kTa 
in  Eq.  2.14  of  this  study).  This  requires  that  the  roots  be  determined  by  numerical 
analysis. 

Brevart  and  Fuller  computed  the  response  of  an  infinite  fluid  filled  elastic 
shell  to  a  radial  impulse  defined  as  a  delta  function  in  both  time  and  space.  The 
axisvmmetric  description  of  this  system  may  be  constructed  using  Eqs.  2.14  and 
3.8  as  follows: 
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where  VW(Q,KZ)  is  the  double  Fourier  transform  of  the  radial  velocity  in 
frequency-wavenumber  space,  ps  is  the  density  of  the  shell  material,  Pj  is  the  den¬ 
sity  of  the  fluid,  cp  is  the  extensional  wave  speed  of  the  shell,  Kr  = 
is  the  non-dimensional  radial  wavenumber  for  the  fluid  cylinder,  and  Pa  is  the 
Fourier  transform  of  the  input  ring  force  function,  P06(t)6(z).  The  authors  ob¬ 
tained  the  time  domain  response  by  numerical  inversion  of  the  three  dimensional 
version  of  Eq.  4.9.  The  inversion  in  wavenumber  space  was  computed  by  sum¬ 
ming  the  residues  of  the  transfer  receptance,  as  previously  described  by  Fuller 
(1983).  This  required  knowledge  of  the  location  of  the  poles  of  the  transfer  mo- 
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bility,  which  occur  at  the  zeros  of  the  characteristic  equation  for  free  vibrations  of 
the  fluid  loaded  cylinder.  The  earlier  paper  by  Fuller  and  Fahy  (1982)  provided 
this  information. 


4.3.3  Elastic  Shell  with  External  Fluid:The  response  of  an  infinite  cylindrical 
shell  with  external  fluid  loading  to  harmonic  excitation  from  a  ring  force  has  been 
studied  by  Photiadis  (1990).  In  contrast  with  the  internal  fluid  problem,  where 
all  radial  impedances  are  purely  reactive  and  energy  flow  is  confined  to  the  axial 
direction  in  the  shell  and  fluid,  the  propagation  of  waves  outward  in  the  external 
fluid  complicates  the  analysis  of  this  case. 

In  terms  of  the  present  study,  the  problem  discussed  by  Photiadis  may  be 
expressed  as 
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in  which  Po{Vt)  is  the  spatial  Fourier  transform  of  the  harmonic  input,  Pc(Q)S(t), 
and  the  Bessel  functions  representing  the  fluid  cylinder  impedance  in  Eq.  4.9  have 
been  replaced  with  the  Hankel  functions  to  represent  the  external  medium,  as 
written  in  Eq.  2.57.  The  other  quantities  in  Eq.  4.10  are  as  defined  above. 

Photiadis  constructs  a  Green’s  function  solution  for  the  radial  velocity  as  a 
function  of  axial  distance  by  inverse  Fourier  transformation  of  the  denominator  of 
Eq.  4.10,  that  is, 


Vw(Q,z)  =  P0(Q)G(z) 
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The  inversion  integral  is  evaluated  by  contour  integration  in  the  complex  Kz 
plane.  The  analysis  is  fairly  recondite  and  will  not  be  pursued  at  length  here. 
Essentially,  Photiadis  evaluates  the  integral  at  the  various  poles  by  introducing 
simplifying  approximations  based  on  physical  reasoning.  The  approximate  analyt¬ 


ical  results  compare  favorably  with  the  exact  values  computed  numerically. 
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CHAPTER  5 
Experimental  Work 

5.1  Introduction 

In  this  chapter  we  discuss  an  experiment  which  was  conducted  to  study  the 
propagation  of  transients  in  a  thin  cylindrical  shell  for  comparison  with  the  the¬ 
oretical  predictions  described  in  Chapter  3.  We  begin  with  a  description  of  the 
experimental  apparatus  and  the  considerations  which  affected  the  design  of  the  ex¬ 
periment.  This  is  followed  by  a  discussion  of  the  measurements  obtained.  Results 
of  the  experimental  work  are  compared  with  the  predicted  values. 

5.2  Design  of  Experiment 
5.2.1  Description  of  Experiment 

The  experimental  apparatus  comprised  a  thin  brass  tube  which  was  freely  sus¬ 
pended  on  very  compliant  supports  (bungee  cords)  at  two  axial  locations.  The 
tube  was  mechanically  driven  by  a  piezoelectric  disk  which  was  bonded  to  one 
end.  The  disk  was  driven  from  an  electronic  function  generator,  thus  various  wave¬ 
forms  of  interest  (such  as  sinusoidal,  triangular,  or  rectangular)  could  be  applied 
to  the  tube.  The  radial  displacement  of  the  tube  was  measured  with  gages  made 
from  poly-vinylidene  fluoride  film  applied  circumferentially  to  the  exterior  wall  of 
the  tube  at  various  locations  of  interest.  Electrical  output  from  the  displacement 
gages  was  amplified,  filtered,  and  displayed  on  a  digital  oscilloscope.  The  oscillo- 
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scope  has  a  transient  capture  capability  and  a  digital  recorder.  Raw  data  obtained 
with  the  digital  oscilloscope  was  saved  on  computer  data  storage  disks  for  off-line 
processing.  A  sketch  of  the  apparatus  is  shown  in  figure  5-1. 

5.2.2  Selection  of  Tube 

The  important  factors  involved  in  the  selection  of  a  tube  for  the  experiment  include 
the  tube  material,  the  length,  the  diameter,  the  wall  thickness,  the  availability,  and 
the  cost.  The  most  significant  factor  in  the  selection  of  the  material  is  the  wave 
speed.  For  a  finite  length  tube,  a  low  wave  speed  maximizes  the  time  in  which 
transients  may  be  observed  before  echoes  occur.  Unfortunately  the  materials  with 
the  lowest  wave  speeds,  such  as  plastics  and  lead,  are  also  lossy.  For  a  given 
material,  the  diameter  will  determine  the  ring  frequency.  This  should  be  chosen 
so  as  to  be  within  the  bandwidth  of  the  available  instrumentation.  The  wall 
thickness  should  be  selected  to  meet  the  criteria  for  a  ’’thin”  shell,  as  indicated  by 
the  analytical  models. 

The  criteria  of  cost  and  availability  limited  the  choices  to  aluminum,  brass, 
copper,  and  steel.  Brass  was  selected  since  it  has  the  lowest  wave  speed  of  the  four 
choices  (a  dilatational  wave  speed  of  4700  meters  per  second,  compared  with  6300 
m/s  for  aluminum,  5000  m/s  for  copper,  and  6100  m/s  for  steel)  and  is  readily 
available.  Standard  brass  tube  is  available  in  12  foot  (3.61  meter)  lengths,  and  this 
would  conveniently  fit  in  the  available  laboratory  space.  A  diameter  of  2  inches  (5 
cm)  was  selected,  which  places  the  ring  frequency  at  approximately  25  kilohertz. 
The  wall  thickness  of  the  experimental  tube  is  0.035  inches  (.89  millimeters),  so 
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the  thickness  to  mean  radius  ratio  is  0.036. 

5.2.3  Displacement  Gages 

In  order  to  observe  wave  propagation  and  dispersion  in  the  experimental  tube,  it 
is  necessary  to  measure  the  components  of  displacement  at  various  axial  positions 
along  the  length  of  the  tube.  Since  axisymmetric  radial  excitation  is  employed,  only 
radial  and  axial  displacements  are  expected.  An  ideal  displacement  gage  would 
have  the  following  attributes:  high  sensitivity,  low  noise,  minimum  interference 
with  the  signal,  and  dimensions  which  are  small  compared  with  the  wavelength 
of  the  stress  waves  to  be  measured.  For  measurement  of  the  radial  displacement, 
a  gage  which  wraps  around  the  circumference  of  the  tube  will  provide  spatial 
integration  of  any  non-axisymmetric  signal  components. 

Graff  (1975)  describes  several  methods  which  have  been  employed  for  the  ex¬ 
perimental  investigation  of  stress  waves.  These  include  the  electrical-resistance 
strain  gage,  capacitance  gages,  inductance  effects,  piezoelectric  devices,  and  photo¬ 
elasticity.  More  recently,  laser  velocimeters  have  been  employed  by  a  number  of 
investigators. 

For  the  present  application,  a  relatively  new  plastic  piezoelectric  material, 
polarized  poly-vinylidene  fluoride  film  (PVDF),  offered  several  advantages.  This 
film  is  produced  by  stretching  extruded  PVDF  film  to  induce  a  phase  transition 
from  the  non-polar  alpha  phase,  to  the  polar  beta  phase.  The  beta  phase  film  is 
then  polarized  by  exposure  to  a  high  electrostatic  field  at  elevated  temperature. 
The  resulting  permanent  alignment  of  molecular  dipoles  in  the  film  produces  the 
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piezoelectric  activity.  Electrodes  are  applied  to  the  film  either  by  sputtering  or 
by  silk-screening  using  conductive  silver  ink.  Electrodes  may  be  produced  in  any 
desired  pattern. 

The  area  density  of  thin  (  52  •  1CT6  meters)  PVDF  film  is  about  93  grams 
per  square  meter,  which  represents  a  negligible  mass  load  when  applied  to  the 
experimental  tube.  The  Young’s  modulus  of  PVDF  is  also  very  low,  2- 109  Newtons 
per  square  meter,  so  that  a  thin  strip  will  have  a  negligible  stiffness  effect  when 
bonded  to  the  tube.  The  film  is  easily  attached  either  with  cement  or  double  sided 
plastic  tape. 

The  electrical  and  mechanical  properties  of  piezoelectric  film  are  categorized  in 
terms  of  a  1,  2,  3  rectangular  coordinate  system.  The  axes  are  oriented  such  that 
1  corresponds  to  length  (the  direction  in  which  the  film  is  stretched),  2  corresponds 
to  the  width,  and  3  to  the  thickness,  which  is  the  direction  of  polarization.  The 
film  is  anisotropic,  so  that  the  piezoelectric  constants  are  greatest  in  the  3,1  and 
3,3  directions,  and  minimum  in  the  3,2  direction.  For  example,  the  g 31  constant  is 
about  1 1  times  greater  than  g^2 ,  which  means  that  the  voltage  produced  by  a  given 
stress  in  the  length  (1)  direction  is  over  20  dB  greater  than  the  voltage  resulting 
from  the  same  stress  in  the  width  (2)  direction.  This  is  advantageous  in  the  current 
application,  as  it  discriminates  against  displacements  which  are  orthogonal  to  the 
desired  displacement. 

Sensitivity  of  PVDF  Gages:  The  PVDF  gages  are  applied  circumferentially  to  the 
exterior  surface  of  the  brass  tube.  The  active  length  of  the  gages  is  155  millimeters, 
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about  5  mm  less  than  the  circumference  of  the  tube.  Since  the  Young’s  modulus 
of  the  tube  is  about  50  times  that  of  the  PVDF  gages,  we  may  assume  that  the 
strain  in  the  gages  is  the  same  as  that  of  the  tube.  The  stress  in  the  gage  is  the 
product  of  the  strain  and  the  Young’s  modulus  of  the  gage,  and  the  electric  field 
Ez  is  gz\  times  the  stress.  The  voltage  developed  is  the  product  of  the  electric  field 
times  the  thickness  of  the  gage.  Thus  the  voltage  V  resulting  from  a  given  strain 
e  is  given  by  V  =  (Ygtggz i)  •  e  ,  where  Yg  is  the  Young’s  modulus  of  the  PVDF 
and  tg  is  the  thickness  of  the  PVDF  film.  Evaluation  of  this  expression  using  the 
values  for  the  circumferential  gages  results  in  V/e  =  22.46  •  103  (volts/unit  strain). 
A  typical  strain  gage  has  a  sensitivity  of  about  5  volts/unit  strain,  thus  the  PVDF 
film  is  about  73  dB  more  sensitive  than  a  strain  gage.  The  radial  displacement  w  is 
eLgj 27T  ,  since  the  length  of  the  gage  Lg  is  approximately  equal  to  the  circumference 
of  the  tube.  Thus  the  radial  sensitivity  of  the  gage  is  w/V  =  1.112  •  10~6  meters 
per  volt,  or  -119  dB  re  1  meter/volt. 

5.2.4  Piezoelectric  Disk  Driver 

Mechanical  input  to  the  experimental  tube  was  generated  with  a  lead  zirconate 
titanate  (PZT)  ceramic  disk  which  was  bonded  to  one  end  of  the  tube  with  epoxy 
cement.  The  disk  is  5.08  cm  in  diameter  and  2.57  mm  thick.  Circular  silver 
electrodes  4.85  cm  in  diameter  are  provided  on  the  flat  faces  of  the  disk.  When  a 
voltage  is  applied  to  the  electrodes,  the  3,1  mode  coupling  in  the  disk  produces  an 
axisymmetric  radial  strain,  which  is  communicated  to  the  brass  tube.  The  static 
radial  displacement  of  a  free  thin  piezoelectric  ceramic  disk  is  given  by:  (Vernitron, 


where  d 3ji  is  the  piezoelectric  strain  constant  for  the  3, 1  mode,  V  is  the  applied 
voltage,  r  is  the  radius  of  the  disk, and  t  is  the  thickness.  For  the  disk  used  in  this 
experiment,  d3jl  =  -109  *  10""12  meters/volt  and  the  ratio  of  radial  displacement 
to  applied  voltage  is  A r/V  =  2.2  •  10~6  millimeters  per  volt. 

The  piezoelectric  disk  will  have  an  infinite  number  number  of  radial  and  thick¬ 
ness  vibrational  modes.  Since  the  disk  is  almost  completely  silvered  on  both  flat 
surfaces,  only  the  n  =  0  (axisymmetric)  radial  modes  will  be  excited.  Since  the 
radius  of  the  disk  is  approximately  10  times  the  thickness,  the  lower  order  ra¬ 
dial  modes  will  be  independent  from  the  thickness  modes.  The  thickness  mode 
frequencies  are  given  by 

=  m  =  1,2,3,... 

(Mason,  1964)  where  cf  is  the  wave  speed  in  the  3,3  direction  at  constant 
electric  displacement  D,  and  m  is  the  mode  number.  Using  the  value  c®  =  4300 
meters/second  (Berlincourt  and  Krueger,  1964)  yields  the  value  f\  =  838  kilohertz 
for  the  first  thickness  mode  of  the  disk. 

The  frequencies  of  the  axisymmetric  radial  modes  are  found  from  the  roots  of 

v  Mv)  =  (1  -vE)  Mv)  (S-1) 

(Meitzler,  O’Bryan,  and  Tiersten,  1973)  where  vE  =  —(s^/ sfi)  =  -315  is  a 
planar  Poisson’s  ratio,  and  the  argument  rj  =  2nf  a/cp.  The  quantity  cp  in  rj  is  a 
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planar  wave  speed, 

sE 

P  [(sfl)2  - 

in  which  p  is  the  density  of  the  ceramic  material,  and  and  sf2  are  the 
elastic  compliance  constants  at  constant  electric  field  E.  For  the  disk  used  in  this 
experiment,  <f  =  3464  meters /second.  Solving  5.1  the  frequencies  of  the  first  seven 
radial  modes  are  found  as  44.61, 117, 186.1,  254.7,  323.1,  391.4  and  459.7  kilohertz, 
respectively. 

Measurements  made  on  the  piezoelectric  disk  prior  to  installation  are  shown 
on  figures  5-2,  5-3,  and  5-4.  Figures  5-2  and  5-3  present  the  magnitude  of  the 
electrical  admittance  of  the  disk  as  a  function  of  frequency,  plotted  on  logarithmic 
axes.  The  peak  of  the  admittance  magnitude  on  figure  5-2  at  44.61  kilohertz 
identifies  the  frequency  of  the  first  radial  mode.  The  admittance  minimum  at  51.33 
kilohertz  occurs  at  the  antiresonant  frequency  for  the  first  radial  mode.  From  these 
frequencies,  the  planar  coupling  factor  was  determined  to  be  kp  —  .55.  (IEEE  1984) 
Radial  modes  two  through  seven  are  easily  identified  on  figure  5-3,  occurring  at 
117.0,  184.4,  250.0,  311.3,  368.6,  and  420.0  kilohertz,  respectively.  The  magnitude 
of  the  ratio  of  maximum  to  minimum  admittance  at  a  specific  mode  decreases  as  the 
mode  number  increases,  thus  it  is  difficult  to  identify  modes  above  the  seventh  on 
the  plot  of  figure  5-3.  The  measured  frequencies  agree  closely  with  those  caclulated 
above,  especially  for  the  first  four  modes.  Note  that  the  difference  between  the 
predicted  and  measured  nodal  frequencies  increased  with  mode  number.  This  may 
be  due  to  the  fact  that  the  diameter  of  the  electrodes  is  23  millimeters  less  than 
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the  diameter  of  the  disk.  This  creates  a  ring  of  inactive  ceramic  at  the  outer  edge 
of  the  disk.  The  mass  of  this  inactive  ring  will  tend  to  lower  the  modal  fequencies, 
becoming  more  significant  at  the  higher  frequencies. 

The  first  thickness  mode  appears  at  808  kilohertz  on  figure  5-3.  The  rapid 
fluctuations  in  admittance  near  and  above  this  resonance  may  be  due  to  modal 
interaction  with  higher  order  radial  modes.  For  example,  the  computed  values  for 
modes  12  to  14  are  800.9,  869.1,  and  937.3  kilohertz,  respectively. 

The  complex  admittance  locus  (susceptance  versus  conductance)  for  the  first 
radial  mode  of  the  disk  is  shown  in  figure  5-4.  The  mechanical  ”Q”,  or  quality 
factor,  of  the  disk  is  determined  from  this  plot  as: 


Q  = 


fr 

(h  ~  h) 


(5.2) 


in  which  fr  is  the  mechanical  resonance  frequency,  the  frequency  at  which  the 
conductance  is  a  mazimum,  and  and  / \  are  the  quadrantal  frequencies,  where 
the  conductance  is  half  of  the  maximum  value.  The  Q  for  the  experimental  disk 
computed  with  equation  5.2  is  816,  which  indicates  that  the  free  disk  has  very  low 
losses. 


5.2.5  Development  of  Dummy  Load 

Initial  measurements  on  the  apparatus  with  a  2  microsecond  rectangular  pulse 
input  applied  to  the  piezoelectric  disk  revealed  the  presence  of  strong  spectral 
components  at  frequencies  near  those  of  the  first  three  radial  modes  of  the  free 
piezoelectric  disk.  Figure  5-5  presents  the  spectrum  of  the  radial  displacement 
at  gage  C-l,  located  approximately  7  cm  from  the  piezoelectric  disk.  Since  it  is 
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desireable  to  have  a  relatively  flat  spectrum  in  the  mechanical  input  to  the  tube, 
various  means  for  achieving  this  end  were  explored. 

The  measurements  on  the  free  piezoelectric  disk,  before  it  was  attached  to  the 
tube,  were  made  with  an  impedance  analyzer  which  applies  a  continuous  harmonic 
excitation  to  the  test  specimen.  However,  driving  the  attached  disk  with  a  con¬ 
tinuous  signal  was  not  feasible  due  to  the  multiple  resonant  modes  which  were 
excited  in  the  tube.  Instead,  the  disk  was  driven  with  a  tone  burst  comprised  of  a 
gated  sinusoidal  voltage  signal  at  the  first  radial  mode,  about  44  kHz.  The  Q  of 
the  system  was  then  estimated  from  the  rise  time  of  the  current  in  the  disk.  With 
this  method,  the  tone  burst  must  be  long  enough  to  encompass  the  initial  and  final 
transients  in  the  system  response. 

Tone  burst  measurements  indicated  that  the  Q  of  the  system  at  the  first  radial 
mode  was  approximately  50.  One  simple  way  to  reduce  the  Q  of  a  system  is  to 
increase  the  damping  by  adding  resistive  losses.  Various  lossy  materials,  such  as 
a  mastic  material,  and  plastic  tape,  were  applied  to  the  free  face  of  the  piezo¬ 
electric  disk,  and  the  free  end  of  the  tube.  The  lowest  value  of  Q  which  could  be 
achieved  with  these  methods  was  about  15,  which  was  still  too  high  for  the  planned 
measurements. 

A  more  sophisticated  approach  to  increasing  the  bandwith  of  piezoelectric 
transducers  is  by  the  use  of  an  electrical  ’’dummy  load.”  With  this  method,  a 
second  piezoelectric  element  is  mechanically  coupled  to  the  driven  device.  The 
electrical  output  of  the  second,  or  ’’dummy,”  element  is  connected  to  a  suitable 
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electrical  network,  which  is  designed  to  maximize  the  bandwidth  of  the  system. 
Application  of  this  method  to  ultrasonic  transducers  has  been  described  by  Gaza- 
let,  Houze,  et  al  (1988).  For  the  present  application,  the  dummy  element  comprises 
a  piezoelectric  disk  which  is  identical  to  the  driven  disk.  The  dummy  disc  is  ce¬ 
mented  to  the  driven  disk,  with  appropriate  means  for  connection  to  an  external 
electrical  network. 

The  first  step  in  the  development  of  the  dummy  load  was  to  prepare  a  test  setup 
comprising  two  piezoelectric  disks  which  were  bonded  together,  but  not  attached 
to  the  tube.  This  enabled  the  dummy  load  system  to  be  optimized  independent  of 
the  additional  load  imposed  when  the  brass  tube  is  driven.  Initial  measurements 
were  made  on  the  two  disks  to  determine  their  free  characteristics,  then  the  disks 
were  cemented  together  to  form  a  prototype  dummy  load  system. 

For  frequencies  in  the  vicinity  of  the  first  radial  mode,  the  electromechanical 
characteristics  of  the  piezoelectric  disk  may  be  represented  by  a  lumped  analogous 
electrical  network,  shown  in  figure  5-6  (a)  (Mason,  1964;  Kinsler,  Frey  et  al,  1982). 
In  figure  5-6  (a)  we  employ  the  impedance  analogy,  in  which  electrical  voltage 
is  analagous  to  force,  and  electrical  current  is  analogous  to  velocity  (Beranek, 
1954).  The  mechanical  elements  (mass,  compliance,  and  loss  resistance)  have  been 
converted  into  equivalent  electrical  circuit  elements  (inductance  Lme ,  capacitance 
Cme,  and  resistance  Rme)  using  the  electromechanical  transformation  factor,  N. 
The  input  electrical  voltage  E  and  current  I  are  also  indicated  on  figure  5-6  (a). 
The  capacitor  Coi  represents  the  ”  blocked”  capacitance  of  the  disk,  that  is,  the 
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capacitance  which  one  would  measure  if  the  disk  was  physically  restrained  from  any 
mechanical  motion  The  mechanical  variables,  force  and  velocity,  are  shown  as  F/N 
and  U  •  N.  All  of  the  elements  shown  in  figure  5-6  (a)  may  be  determined  from  the 
results  of  two  terminal  measurements  of  a  ’’free”  (that  is,  mechanically  unloaded) 
disk.  Since  air  has  a  very  low  specific  acoustic  impedance,  the  mechanical  terminals 
of  figure  5-6  (a)  are  effectively  ’’short  circuited”  for  air  measurements,  resulting 
in  the  equivalent  cirtcuit  of  figure  5-6  (b).  First,  the  capacitance  of  the  disk  is 
measured  at  a  frequency  far  below  resonance.  This  ’’free”  capacitance  is  equal 
to  the  parallel  combination  of  Co  and  Cme,  which  in  turn  may  now  be  computed 
using  the  measured  capacitance  and  the  planar  coupling  factor,  kP.  The  remaining 
elements,  Lme  and  Rrne  are  found  using  the  complex  admittance  locus,  such  as 
shown  in  figure  5-4.  The  mechanical  resonance  frequency  fr  is  the  frequency  at 
maximum  conductance,  Gmax.  The  inductance  Lme  =  1/Cme(2irfr)2,  and  the 
resistance  Rme  —  1/Gmax.  These  quantities  were  determined  from  measurements 
of  several  disks  to  implement  the  next  stage  in  the  dummy  load  design. 

Figure  5-6  (c)  represents  the  equivalent  circuit  of  the  dummy  load  system. 
The  driven  disk  is  represented  by  the  elements  with  subscript  1,  and  the  dummy 
disk  by  elements  with  subscript  2.  The  blocks  Zme\  and  Zme 2  comprise  the  se¬ 
ries  mechanical  branch  of  each  disk.  The  dummy  disk  is  connected  to  a  circuit 
comprising  a  resistor  Rd  in  parallel  with  an  inductor  La-  When  the  dummy  load 
system  is  connected  to  a  mechanical  load,  such  as  the  experimental  brass  tube,  a 
load  impedance  Zie  would  appear  at  terminals  a-a’  in  figure  5-6  (c). 
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A  mathematical  model  was  made  of  the  equivalent  circuit  of  figure  5-6  (c)  using 
average  values  for  the  circuit  elements  from  the  measurements  described  above. 
The  model  enabled  prediction  of  the  response  of  the  dummy  load  system  as  function 
of  the  electrical  dummy  load  resistor  and  inductor.  By  extensive  exercise  of  the 
model  it  was  determined  that  a  Q  of  2.1  could  be  achieved  with  the  proper  values 
of  Rd  and  Ld.  Using  these  results,  input  admittance  measurements  were  made  on 
the  prototype  dummy  load  system  with  various  combinations  of  Rd  and  Ld-  The 
lowest  Q  obtained  experimentally  on  the  prototype  system  was  4,  determined  from 
analysis  of  the  complex  admittance  locus. 

The  next  effort  was  focused  on  the  application  of  the  dummy  load  to  the 
experimental  apparatus.  A  dummy  load  disk  was  bonded  to  the  primary  driving 
disk  with  epoxy  cement.  A  thin  (0.2  millimeter)  circular  disk  of  expanded  metal 
was  placed  between  the  piezoelectric  disks  to  facilitate  the  electrical  connection 
to  the  common  silver  interfaces  of  the  disks.  The  dummy  load  disk  was  then 
connected  to  the  electrical  dummy  load  which  had  been  developed  previously  using 
the  prototype  dummy  load  system. 

Addition  of  the  dummy  load  produced  a  dramatic  improvement  in  the  opera¬ 
tion  of  the  experimental  system.  This  is  illustrated  in  figure  5-7,  which  compares 
the  system  response  with  the  dummy  load  ’’off,”  that  is,  electically  disconnected, 
with  that  when  the  electrical  dummy  load  was  connected  to  the  dummy  load  disk. 
The  primary  disk  was  driven  with  a  3.3  microsecond,  30  volt  rectangular  pulse  for 
these  measurements.  The  output  was  taken  as  the  voltage  developed  across  the 
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dummy  load  disk,  which  is  proportional  to  the  radial  displacement. 

Figure  5-7  (a)  displays  the  dummy  disk  voltage  with  the  dummy  load  ’’off.” 
The  system  ’’rings”  for  about  800  microseconds  -with  this  input.  The  spectrum  of 
this  waveform  is  shown  in  figure  5-7  (c),  in  which  a  logarithmic  scale  is  used  on 
the  ordinate.  The  major  peak  in  the  response  occurs  at  47.8  kilohertz,  and  the  Q 
as  determined  from  the  frequencies  at  which  the  response  has  decreased  by  3  dB 
is  about  28. 

When  the  dummy  load  is  connected  (figure  5-7  (b)),  the  response  decays  to 
zero  in  less  than  100  microseconds.  The  spectrum,  in  figure  5-7  (d),  peaks  at  43 
kilohertz.  The  Q  has  dropped  to  2.7,  a  ten  to  one  improvement  over  the  first 
case.  Figure  5-7  (d)  also  shows  that  the  spectrum  from  100  to  200  kilohertz  has 
been  elevated  about  20  dB,  indicating  that  the  dummy  load  also  improves  the  high 
frequency  response  of  the  system. 

5.2.6  Determination  of  Wave  Speed 

In  order  to  insure  an  accurate  comparison  of  the  results  of  experimental  measure¬ 
ments  with  those  predicted  by  theory,  it  is  necessary  to  know  the  exact  value  of  the 
wave  speed  in  the  experimental  brass  tube.  Five  independent  methods  were  inves¬ 
tigated  for  this  purpose.  First,  a  review  was  made  of  published  values  for  brass. 
The  second  approach  was  a  calculation  from  the  material  properties  published  for 
the  type  of  brass  used.  Then  three  experimental  methods  were  considered:  (1) 
Time  of  Flight  measurement;  (2)  Ring  Mode  Resonance;  and  (3)  Extrapolation 
from  free  longitudinal  modal  resonances.  The  most  accurate  result  was  attained 
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with  the  last  method. 

Review  of  Published  Data:  The  results  of  a  review  of  three  texts  Fahy  (1985  ), 
Elmore  and  Heald  (1985  ),  and  Kinsler  and  Frey,  et.  al.  (1982  ))  and  tables 
from  a  commercial  source  (Metrotek  (1982))  are  presented  in  Table  5-1.  The  entry 
labeled  ”  Cartridge  Brass,”  which  is  the  type  of  brass  for  the  experimental  tube,  was 
computed  from  values  of  Young’s  modulus  and  density  from  Hoyt  (1954).  The  table 
lists  the  source,  the  value  of  Young’s  Modulus  E,  in  Newtons  per  square  meter;  the 
density  p ,  in  kilograms  per  cubic  meter;  Poisson’s  ratio  v,  the  longitudinal  wave 
speed  (also  known  as  the  ’’bulk  wave  speed”  and  the  ’’dilatational  wave  speed”) 
cd:  the  ’’thin  plate”  (or  ’’extensional”)  wave  speed  cp;  the  ’’thin  bar”  wave  speed 
Cj,;  the  shear  wave  speed  cs;  (all  wave  speeds  are  in  units  of  meters  per  second);  the 
shear  modulus  of  elasticity  G  in  Newtons  per  square  meter  (this  is  also  the  Lame 
coefficient  /i);  the  bulk  modulus  B  in  Newtons  per  square  meter;  and  the  Lame/ 
coefficient  A,  Newtons  per  square  meter. 

None  of  the  sources  reviewed  published  all  of  the  properties  listed  in  each  row 
of  Table  5-1.  Any  missing  values  were  calculated  from  those  published,  using  one 
or  more  of  the  equations  listed  below. 

Dilatational  (longitudinal)  wave  speed: 

I J _  E{l-u) 

Cd  ~  \  p  ~  \  p{l  +  v)  (l  -  2v) 


Extensional  (thin  plate)  wave  speed: 
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Thin  bar  wave  speed: 


c‘  ”  V  p  ~  M 


fi  (3A  +  2  /jl) 
p(A  +  /i) 


Shear  wave  speed: 


cs  — 


E  Hx 


2p  (1  +  v)  V  P 


Bulk  modulus  of  elasticity: 


f3  —  \  +  -fi 


Young’s  modulus  of  elasticity: 


E  —  2/i  (1  1/)  — 


(3A  2/x) 

A  +  /i 


Poisson’s  ratio: 

_  A 
2  (A  +  /i) 

The  mean  and  standard  deviation  of  the  published  data  are  shown  in  the  table. 
The  deviation  expressed  as  a  percentage  of  the  mean  is  also  included.  The  large 
variation  in  the  published  data  on  brass  may  be  due  to  the  fact  that  brass  is  an 
alloy  of  copper  and  zinc,  with  trace  quantities  of  lead  and  iron.  Hoyt  (1954)  lists 
over  13  different  types  of  brass.  The  copper  content  varies  from  60%  to  80%, 
and  the  lead  from  .07%  to  2%.  The  iron  content  ranges  from  .05%  to  .15%.  The 
balance  of  the  alloy  is  zinc.  The  mean  value  of  Young’s  modulus  for  these  alloys 
is  10.65  X  1010  Newtons/square  meter,  with  a  standard  deviation  of  5.76%,  and 
the  mean  value  of  density  is  8.47  x  103  kg/cubic  meter,  with  a  standard  deviation 
of  1.16%.  These  values  compare  closely  with  those  in  Table  5-1.  Unfortunately, 
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Hoyt  does  not  list  the  values  of  Poisson’s  ratio.  The  composition  of  the  brass  tube 
used  in  the  experiment  is:  copper,  68.5%  minimum,  71.5%  maximum;  lead,  0.07%, 
iron,  0.05%,  and  the  remainder  is  zinc.  This  composition  fits  the  description  of 
cartridge  brass  in  Hoyt. 

Time  of  Flight  Measurement:  In  the  time  of  flight  method,  the  wave  speed  is 
measured  by  determination  of  the  time  required  for  a  transient  wave,  such  as 
an  impulse  or  short  rectangular  pulse,  to  traverse  a  known  distance.  For  a  non- 
dispersive  system  (such  as  a  string,  or  long  thin  bar),  this  method  yields  accurate 
results.  With  a  dispersive  system,  such  as  the  experimental  tube,  distortion  of  the 
input  pulse  with  distance  makes  it  difficult  to  identify  the  precise  time  of  arrival 
of  the  signal.  Nine  measurements  on  the  experimental  tube  gave  values  of  cv 
ranging  from  3696  to  4507  meters  per  second.  The  average  value  was  3959  meters 
per  second,  with  a  standard  deviation  of  244  meters  per  second.  These  results 
indicated  the  need  for  a  more  precise  method  for  the  determination  of  wave  speed. 

Ring  Mode  Resonance:  Another  potential  method  for  determining  cv  is  based 
on  harmonic  excitation  of  the  so-called  "ring”  frequency  of  the  tube.  With  an 
axisymmetric  (n  =  0  mode)  excitation,  a  strong  resonance  should  occur  at  the 
"ring”  frequency,  where  the  axial  wavelength  is  equal  to  the  mean  circumference 
of  the  tube.  The  extensional  wave  speed  may  then  be  computed  from  the  relation 
cp  =  27 vfra,  where  fr  is  the  measured  ring  frequency,  and  a  is  the  mean  radius 
of  the  tube.  This  measurement  was  attempted  by  applying  continuous  sinusoidal 
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excitation  to  the  piezoelectric  driver  disk  and  monitoring  the  radial  displacement 
of  the  tube  on  the  circumferential  PVDF  gages.  Using  the  estimated  value  of 
cp  =  3872  meters/sec  from  Table  W-l,  the  ring  frequency  should  occur  at  24.69 
kilohertz. 

The  problem  with  this  approach  is  that  since  a  finite  length  tube  is  used,  many 
longitudinal  modes  are  excited  when  a  continuous  harmonic  excitation  is  applied. 
Modes  were  found  at  24.431  kHz,  24.652  kHz,  and  25.811  kHz.  The  dilemma  in 
this  approach  is  in  separating  the  true  ring  frequency  from  the  longitudinal  modal 
frequencies.  This  difficulty,  in  fact,  suggested  that  the  low  frequency  longitudinal 
modes  might  be  used  to  more  accurately  asses  the  wave  speed  of  the  tube. 

Longitudinal  Modal  Resonance:  This  method  is  based  on  an  analysis  of  a  finite 
length  tube  driven  at  frequencies  far  below  the  ’’ring”  frequency,  so  that  bending 
stresses  may  be  neglected.  The  experimentally  determined  modal  frequencies, 
which  depend  on  the  tube  dimensions  and  the  wave  speed,  are  compared  with 
those  predicted  by  the  model  developed  in  Chapter  3.  The  value  of  wave  speed  in 
the  model  is  adjusted  until  close  agreement  between  the  predicted  and  measured 
modal  frequencies  is  achieved. 

The  modal  resonance  frequencies  were  measured  on  the  experimental  tube  for 
the  first  eight  longitudinal  modes.  The  piezoelectric  disk  was  driven  with  a  contin¬ 
uous  sinusoidal  waveform  from  a  frequency  synthesizer  which  could  be  controlled  in 
increments  of  0.1  Hertz.  The  radial  and  longitudinal  displacements  were  monitored 
on  four  PVDF  gages.  The  approximate  location  of  the  modal  frequencies  was  de- 
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termined  with  equation  3.25  using  the  value  of  cv  —  3872  meters/second  calculated 
for  cartridge  brass.  For  a  particular  mode,  the  frequency  synthesizer  was  set  to 
the  calculated  frequency  for  that  mode.  The  frequency  was  then  adjusted  until  a 
maximum  (peak  displacement)  was  observed  on  one  of  the  four  gages.  (Due  to  the 
sinusoidal  nature  of  the  mode  shapes,  a  particular  gage  may  be  at  a  displacement 
node  for  some  of  the  modes.) 

The  measured  modal  frequencies  are  shown  in  Table  5-2,  which  also  includes 
the  modal  frequencies  computed  from  the  model.  Table  5-2  also  includes  for  com¬ 
parison  a  list  of  modal  frequencies  computed  for  a  free-free  bar  of  the  same  length 
and  material  as  the  experimental  tube.  Such  a  bar  would  be  non-dispersive,  with 
modal  frequencies  given  by 


j" ndm 


m  •  CgyJ  1  —  z/2 


where  Cpy/1  -  v2  =  cb  =  ^  is  the  thin  bar  wave  speed. 

The  best  estimate  for  the  extensional  wave  speed  was  determined  by  an  itera¬ 
tive  process  in  which  a  range  of  values  were  assumed.  For  each  assumed  value,  the 
modal  frequencies  were  computed  for  the  membrane  model,  and  also  for  the  free- 
free  bar  model.  Then  the  difference  between  the  measured  and  calculated  modal 
frequencies  was  obtained  for  the  first  eight  modes.  A  typical  plot  of  the  difference 
frequency  as  a  function  of  mode  number  is  shown  in  figure  5-8.  A  positive  value  for 
the  difference  indicates  that  the  calculated  modal  frequency  exceeds  the  measured 
modal  frequency  for  that  mode. 

Then  the  mean  and  standard  deviation  was  computed  for  each  set  of  differ- 
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ence  frequencies.  The  standard  deviation  is  plotted  in  figure  5-9  as  a  function  of 
the  assumed  values  of  extensional  wave  speed.  Figure  5-9  indicates  a  minimum 
standard  deviation  of  2  Hertz  for  the  membrane  model  at  a  wave  speed  of  3939  me¬ 
ters/second.  Thus  this  value  was  chosen  as  being  the  most  likely  from  a  statistical 
viewpoint. 

Figure  5-10  presents  the  difference  between  the  computed  and  measured  modal 
frequencies  expressed  as  a  percentage  of  the  measured  frequency,  using  the  final 
value  of  3939  meters /second.  The  maximum  deviation  is  in  the  order  of  0.2%. 

5.3  Transient  Measurements 

As  indicated  earlier,  transient  waves  were  applied  to  the  experimental  tube  by 
exciting  the  piezoelectric  disk  with  appropriate  electrical  waveforms,  such  as  a  short 
rectangular  pulse.  Circumferential  PVDF  gages  were  placed  at  axial  locations  3.3 
(Gage  C-3),  7.3  (Gage  C-l),  19.6  (Gage  C-2),  and  151  (Gage  C-2)  centimeters. 
Gage  outputs  were  applied  directly  to  adjustable  gain,  high  impedance  instru¬ 
mentation  amplifiers.  Filtered  outputs  from  the  amplifiers  were  observed  with  a 
digital  oscilloscope  which  had  an  8  bit,  20  megahertz  analog  to  digital  converter. 
Captured  signals  were  stored  in  a  2  by  4096  memory,  from  which  they  could  be 
output  to  a  computer  via  an  RS-232  interface.  The  digital  oscilloscope  placed  two 
constaints  on  the  quality  and  quantity  of  data  which  could  be  measured.  The  8  bit 
A/D  converter  limits  the  resolution,  and  the  memory  length  of  4096  points  limits 
the  length  of  record  which  may  be  stored.A  variety  of  electronic  signal  generators 
were  used  to  generate  input  to  the  system. 
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The  finite  length  of  the  tube,  approximately  361  centimeters,  meant  that  mul¬ 
tiple  echos  would  be  present  in  long  records.  The  fastest  wave  components  travel 
at  0.3939  cm/microsecond,  so  the  first  echo  at  a  particular  gage  may  be  expected 
at  te  =  (722  -  Lg)/(. 3939)  microseconds,  where  Lg  is  the  distance  from  the  input 
disk  to  the  gage.  Thus  for  gage  C-2,  at  151  cm,  the  echo  appears  about  1450 
microseconds  after  the  first  arrival  at  gage  C-2. 

A  large  number  of  measurements  were  made  over  the  course  of  the  study. 
One  such  series  will  be  described  to  illustrate  the  results  obtained.  In  this  set  of 
tests,  the  input  disk  was  driven  with  a  rectangular  pulse  about  3.3  microseconds 
in  length. 

The  waveforms  observed  at  the  four  gages  are  reproduced  in  figure  5-11.  The 
effects  of  dispersion  are  evident  in  this  figure,  as  the  length  of  the  transient  increases 
at  greater  axial  distances.  The  structure  of  the  waveforms  is  similar  to  those 
which  were  computed,  especially  for  the  first  two  signals,  exhibiting  a  burst  of 
high  frequency  energy  followed  by  a  long  tail  of  long  period  waves.  Echos  are 
evident  on  gages  C-4,  at  19.6  cm,  and  C-2,  at  151  cm.  The  echo  on  gage  C-2  is 
compressed  in  amplitude  and  extended  in  time  compared  to  the  direct  arrival.  The 
time  scale  for  figure  5-11  is  in  milliseconds. 

The  frequency  spectra  for  the  waveforms  in  figure  5-11  are  shown  in  figure 
5-12,  in  which  the  magnitude  of  the  FFT  (normalized  to  the  level  of  the  largest 
component)  is  plotted  as  a  function  of  frequency  in  kilohertz.  The  spectrum  for 
gage  C-3  at  3.3  cm  is  similar  to  that  discussed  earlier  for  the  dummy  load  disk. 
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At  gage  C-l  (7.3  cm),  the  spectrum  is  more  complex,  as  the  transient  moves  down 
the  tube.  The  echos  arriving  at  gages  C-4  (19.6  cm)  and  C-2  (151  cm)  produce 
rapid  fluctuations  in  the  spectra,  figures  5-12  (c)  and  5-12  (d). 

The  effects  of  dispersion  are  clearly  evident  in  the  spectograms,  figures  5-13  to 
5-16.  At  gage  C-3  (fig  5-13))  we  see  the  lag  in  the  wave  components  near  the  ring 
frequency  beginning  to  form  as  a  little  ”nub.”  By  the  time  the  wave  reaches  gage 
C-l  (figure  5-14)  a  well  developed  ’’tail”  of  low  frequency  energy  is  evident. 

Echos  complicate  the  spectograms  at  gages  C-4  and  C-2.  Figure  5-15  clearly 
shows  dispersion  on  the  first  arrival,  and  the  presence  of  the  echo  at  2  milliseconds. 
The  lack  of  a  ”  tail”  on  the  echo  is  due  to  the  limited  length  of  the  time  series  (see 
fig  5-11  (c))  which  has  truncated  the  echo.The  echo  and  the  first  arrival  are  closer 
in  time  for  gage  C-2  (figure  5-11  (d))  and  thus  more  of  the  echo  energy  appears 
in  the  spectogram,  figure  5-16..  The  ’’boomerang”  shape  of  the  joint  distribution 
is  clearly  evident  on  this  figure.  The  echo  also  displays  the  increased  spreading 
in  both  frequency  and  time  as  a  consequence  of  the  transient  travelling  a  greater 
distance  than  the  first  arrival. 

Further  insight  is  provided  by  the  three  dimensional  contour  plots  of  figures  5- 
17  to  5-20.  In  particular,  the  magnitude  of  the  components  of  the  joint  distribution 
is  more  apparent  in  the  3  D  contours.  The  peaks  in  the  frequency  domain  correlate 
well  with  the  maximum  in  the  spectra  of  figure  5-12.  The  growth  of  the  ring 
frequency  ’’tail”  is  evident  in  figures  5-17  and  5-18.  The  stucture  of  the  first 
arrival  and  the  echo  at  gagaes  C-4  and  C-2  is  evident  in  figures  5-19  and  5-20. 
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5.3.1  Comparison  of  Measured  and  Computed  Results 

The  theoretical  results  presented  in  Chapter  3  were  obtained  with  an  idealized 
band-limited  frequency  spectrum.  In  order  to  apply  the  analysis  to  the  experi¬ 
mental  situation,  it  is  first  necessary  to  simulate  the  spectrum  of  the  test  input. 
A  first-order  approximation  to  the  experimental  input  spectrum  shown  in  figure 
5-7  (d)  is  given  by  the  Gaussian  frequency  distribution  described  in  equation  3.53 
using  a  cutoff  frequency  Fc  of  45  kilohertz  and  a  standard  deviation  of  1.7.  The 
predicted  waveform  of  the  radial  displacement  at  gage  C-l,  located  at  7.3  cm., 
is  shown  in  figure  5.21.  The  measured  displacement  is  presented  in  figure  5-11 
(b).  The  two  waveforms  are  similar  in  general  aspects.  The  spectogram  for  the 
computed  waveform  shown  in  figure  5-22  may  be  compared  with  the  experimental 
spectogram  at  gage  C-l,  figure  5-14.  The  tip  of  the  trailing  tail  occurs  at  23.4 
kHz  in  the  caculated  spectrogram,  and  at  23.75  kHz  on  the  experimental  plot. 
Both  figures  indicated  a  concentration  of  energy  at  38  kHz,  and  around  55  kHz. 
The  maximum  frequency  component  occurs  at  89  kHz  in  figure  5-22,  versus  83.6 
kHz  for  the  experiment. 

The  major  difference  is  in  the  overall  width  in  time.  The  computed  spectogram 
is  about  0.4  milloseconds  in  width,  whereas  the  experimental  spectogram  shows 
a  1  milllisecond  duration.  The  discrepancy  is  probably  due  to  the  approximate 
nature  of  the  Gaussian  simulation  to  the  experimental  spectrum. 
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CHAPTER  6 
Conclusions 

6.1  Transfer  Matrix  Method 

The  transfer  matrix  method  is  a  simple  and  straightforward  approach  to  the 
solution  of  wave  propagation  problems  for  concentric  cylindrical  systems  with  fluid 
and  elastic  components.  Models  were  developed  for  a  fluid  cylinder,  a  fluid  tube, 
a  thin  elastic  cylinder,  and  the  external  fluid.  These  models  are  easily  combined 
to  represent  a  complex  system  comprised  of  any  number  of  components. 

Practical  applications  of  this  method  include:  (1)  Analysis  of  cylindrical  sonar 
transducers,  (2)  Study  of  the  effects  of  turbulent  boundary  layer  noise  on  the 
performance  of  towed  arrays,  (3)  Analysis  of  sound  transmission  and  noise  reduc¬ 
tion  in  aircraft,  submarines,  piping  systems,  and  other  structures  which  may  be 
approximated  by  elastic  cylinders. 

The  fluid  models  developed  in  this  work  were  studied  to  examine  the  prop¬ 
agation  of  free  waves  in  the  fluid  cylinder,  the  fluid  tube,  and  the  external  fluid 
for  fixed  and  free  boundary  conditions.  These  limiting  values  are  useful  in  under¬ 
standing  the  effects  of  fluid  density,  sound  speed,  and  the  characteristics  of  external 
excitation  (frequency  and  axial  wavenumber)  in  multicomponent  systems. 

Several  examples  of  composite  systems  are  described  in  Chapter  4,  which  be¬ 
gins  with  a  discussion  of  the  parameters  for  two  fluid  problems  and  the  conditions 
under  which  wave  propagation  may  occur.  A  detailed  study  of  a  fluid  tube  in  an 


137 


external  fluid  illustrates  the  application  of  the  transfer  matrix  method,  and  calcu¬ 
lations  and  plots  of  the  radial  spectral  impedance  at  the  inner  surface  of  the  fluid 
tube  are  shown  for  various  combinations  of  fluid  tube  thickness  and  the  charac¬ 
teristics  of  the  two  fluids.  A  practical  application  in  which  a  fluid  tube  is  used  to 
match  the  impedance  of  a  source  to  that  of  the  external  fluid  is  described. 

The  special  case  of  a  system  which  includes  a  fluid  cylinder  is  described.  It  is 
shown  that  the  fluid  cylinder  has  a  major  influence  on  the  radial  spectral  impedane 
of  the  system. 

Two  examples  of  systems  including  elastic  shells  are  described  briefly. 

6.2  Elastic  Media 

An  extensive  study  of  the  response  of  a  semi-infinite  thin  elastic  cylinder  to 
axisvmmetric  harmonic  and  transient  excitation  is  presented.  Significant  results 
include  the  transfer  matrix  model  for  a  thin  shell,  an  accurate  method  for  the 
determination  of  extensional  wave  speed  in  a  harmonically  driven  finite  length 
shell,  the  effect  of  shell  thickness  on  the  roots  of  the  dispersion  equation,  and  the 
dispersion  of  transient  waves  in  a  thin  shell.  Joint  time-frequency  distributions 
are  shown  to  be  of  value  in  understanding  the  propagation  of  transient  waves  in 
dispersive  media. 

6.2.1  Thin  Shell  Model 

The  transfer  matrix  model  of  a  thin  shell  is  developed  for  use  in  modeling  mul¬ 
ticomponent  cylindrical  systems.  This  model  is  valid  for  shells  which  satisfy  the 
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assumptions  made  in  the  formulation  of  the  Donnell  shell  equations.The  radial 
spectral  impedance  of  a  thin  shell  is  shown  to  be  a  function  of  the  excitation 
frequency,  the  axial  wavenumber,  the  shell  thickness  and  radius,  and  the  den¬ 
sity,  elastic  modulus  and  Poisson’s  ratio  of  the  shell  material.  The  radial  spectral 
impedance,  that  is,  the  ratio  of  the  differential  pressure  across  the  shell  to  the  shell 
radial  velocity,  is  the  principal  term  in  the  two-port  transfer  matrix  for  the  thin 
shell. 

6.2.2  Determination  of  Wave  Speed 

A  very  accurate  method  was  developed  for  the  determination  of  the  extensional 
wave  speed  in  a  thin  elastic  shell.  In  this  approach,  the  longitudinal  resonant 
frequencies  are  measured  for  the  lower  order  modes  of  the  finite  length  shell.  The 
wave  speed  is  computed  from  these  frequencies  using  a  model  based  on  membrane 
theory,  which  is  accurate  for  modal  frequencies  which  are  below  the  ring  frequency. 
The  extensional  wave  speed  for  the  experimental  shell  was  3939  meters/second, 
with  a  maximum  estimated  error  in  the  order  of  0.2%. 

6.2.3  Effect  of  Shell  Thickness  on  the  Dispersion  Equation 

The  dispersion  equation  for  a  thin  shell  includes  a  bending  term,  j3,  which  is 
proportional  to  the  ratio  of  shell  thickness  to  radius.  When  j3  is  zero,  the  shell  is 
characterized  by  membrane  theory  and  there  are  two  real  roots  in  the  solution  of 
the  dispersion  equation.  These  roots  define  a  lower  transition  frequency,  at  which 
free  wave  propagation  of  flexural  waves  in  a  membrane  is  cut  off,  and  an  upper 
transition  frequency,  at  which  longitudinal  wave  propagation  cuts  on.  The  upper 
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transition  frequency  is  shown  to  be  the  traditional  ’’ring”  frequency  of  the  shell. 
The  separation  of  the  transition  frequencies  depends  solely  on  the  Poisson’s  ratio 
for  the  material. 

When  f3  is  finite,  the  dispersion  equation  is  of  sixth  order  in  terms  of  the  axial 
wavenumber.  The  bending  term  causes  the  dispersion  curve  of  the  flexural  branch 
to  dip  to  a  minimum  and  then  to  increase  indefinitely  (a  more  sophisticated  theory 
which  includes  rotatory  inertia  and  shear  effects  is  required  to  correctly  predict  the 
flattening  of  the  dispersion  curve  at  very  high  frequencies).  The  bending  effects  are 
responsible  for  the  propagation  of  very  slow  waves  for  frequencies  in  the  transition 
zone. 

6.2.4  Transfer  Functions  in  Frequency  and  Spatial  Domains 

Using  Fourier  and  Laplace  transform  methods,  two  dimensional  transfer  functions 
for  the  radial  and  axial  displacements  were  constructed  with  frequency  and  axial 
distance  as  arguments.  These  transfer  functions  are  studied  to  provide  understand¬ 
ing  of  the  propagation  of  waves  resulting  from  arbitrary  axysimmetric  inputs  to 
the  shell.  The  near  field  and  far  field  behavior  of  these  shell  displacements  are 
inferred  from  two  and  three  dimensional  representations  of  the  transfer  functions. 

6.2.5  Time  and  Frequency  Domain  Response  of  the  Thin  Shell 

The  radial  and  axial  displacement  responses  of  the  shell  as  functions  of  time  and 
frequency  at  arbitrary  axial  positions  are  evaluated  by  inverse  Fourier  transfor¬ 
mation  of  the  transfer  functions  on  the  frequency  variable.  A  weighting  function 
technique  is  introduced  to  band  limit  the  transfer  functions  in  frequency.  This  is 


140 


shown  to  be  required  to  prevent  aliasing  when  using  discrete  versions  of  the  inverse 
Fourier  transform,  and  also  due  to  limitations  on  the  bandwidth  of  the  transfer 
functions  which  are  imposed  by  the  simplifications  inherent  in  the  formulation  of 
the  Donnell  thin  shell  theory.  The  weighting  function  may  also  be  tailored  to  rep¬ 
resent  the  spectra  of  various  special  temporal  inputs,  such  as  discrete  frequency 
components. 

Plots  of  the  radial  and  axial  displacements  as  functions  of  time  clearly  show 
the  dispersive  nature  of  the  waves  propagating  in  a  thin  shell.  Typically,  the  high 
frequency  components  arrive  first,  followed  by  a  long  "tail”  of  slower  compononets 
at  frequencies  in  the  transition  zone.  The  length  of  the  transient  increases  at 
greater  axial  distances,  due  to  the  spatial  separation  of  the  fast  and  slow  compo¬ 
nents  of  the  transient  wave.  A  beating  effect  is  evident  in  the  displacements  due  to 
constructive  and  destructive  intereference  between  nearby  frequency  components 
in  the  transition  zone. 

The  temporal  variation  of  the  spectral  components  of  the  transient  displace¬ 
ments  are  clearly  shown  by  two  and  three  dimensional  contour  plots  of  the  joint 
time-frequency  distribution  obtained  by  application  of  the  spectogram  to  the  time 
waveform.  The  spectograms  are  shaped  like  a  boomerang  in  the  t  —  f  plane,  with  a 
long  tail  parallel  to  the  time  axis  representing  the  slow,  low  frequency  components. 

6.3  Experimental  Work 

Chapter  5  describes  an  apparatus  constructed  to  generate  and  observe  tran¬ 
sient  wave  propagation  in  a  thin  cylindrical  shell.  Gages  made  from  polarized 
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poly-vinylidene  fluoride  file  (PVDF)  were  found  to  be  an  excellent  means  for  the 
measurement  of  radial  displacement.  These  gages  have  high  sensitivity  and  do  not 
provide  a  significant  mechanical  load  on  the  metallic  cylinder. 

Mechanical  input  to  the  shell  was  provided  by  an  electrically  driven  piezoelec¬ 
tric  ceramic  disk.  An  electro-mechanical  driver  has  several  advantages  over  strictly 
mechanical  devices  (such  as  blows  from  an  instumented  hammer,  explosive  charges, 
or  the  impact  from  projectiles).  Precise,  repetitive  inputs  of  arbitrary  waveform 
are  easily  obtained  from  an  electronic  function  generator.  This  approach  enables  a 
large  volume  of  experimental  data  to  be  obtained  in  a  relatively  short  time  span. 
The  ability  to  accurately  repeat  the  test  stimuli  enhances  the  reliability  of  the 
results  as  compared  to  the  less  controllable  mechanical  sources. 

The  initial  test  system  had  very  low  losses,  which  limited  the  available  band¬ 
width.  This  was  corrected  by  use  of  an  electomechanical  dummy  load  technique,  in 
which  a  second,  identical  piezoelectric  ceramic  disk  was  bonded  to  the  electrically 
driven  disk.  The  second  disk  was  connected  to  an  electrical  network  which  pro¬ 
vided  overall  damping  to  the  system.  The  result  was  a  ten  fold  increase  in  system 
bandwidth. 

Extensive  measurements  of  transient  wave  propagation  were  made  with  the 
experimental  apparatus.  The  results  are  in  good  agreement  with  those  predicted 
by  the  theory.  The  joint  time-frequency  distribution  analysis  was  shown  to  be  an 
effective  tool  in  understanding  wave  propagation  in  dispersive  systems. 
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Classification  of  Cylindrical  Models 
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Fig.  1-7  Geometry  for  Synthesis  of  Cylindrical  Wave 


Fig.  2-2  Modified  Quotient  of  Bessel  Functions  vs. 
non-dimensional  radial  wavenumber,  Kra 
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Fig.  2-4  Normalized  spectral  impedance  of  fluid  cylinder 
as  a  function  of  non-dimensional  wavenumber  Ka,  for  Kz  =  0 
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(b)  Group  Velocity 


Fig.  2-5  Normalized  phase  and  group  velocities  vs  non-dimensional  frequency 
for  a  fluid  cylinder  with  free  radial  boundary,  modes  1,  2,  and  3 
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Fig.  2-6.  Normalized  pressure  vs  non-dimensional  radius 
for  a  fluid  cylinder  with  free  radial  boundary. 
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(a)  Phase  Velocity 


Fig.  2-7.  Normalized  phase  and  group  velocities  vs  non-dimensional  frequency 
for  a  fluid  cylinder  with  fixed  radial  boundary,  modes  1,  2,  and  3. 


Fig.  2-1 1 ,  Surface  Plot  of  Parameter  A 


Fig.  2-12,  Surface  Plot  of  Parameter  B 
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Fig.  2-13,  Surface  Plot  of  Parameter  C 


Fig.  2-14,  Surface  Plot  of  Parameter  D 
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Fig.  2-15  Percent  error  in  the  asymptotic  expression  for  the 
roots  of  the  free-free  mode  in  a  fluid  tube  vs  mode  number. 


Fig.  2-16  Percent  error  in  the  asymptotic  expression  for  the 
roots  of  the  fixed-free  modes  in  a  fluid  tube  vs  mode  number. 


Fig.  2-19  Characteristic  equations  for  free-free  and  fixed-fixed 
modes  of  a  fluid  tube  for  radii  ratios  s  =  1 .5,  2,  and  3 
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Fig.  2-20  Characteristic  equations  for  free-fixed  and  fixed-free 
modes  of  a  fluid  tube  for  radii  ratios  s  =  1.5, 2,  and  3 
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(b)  Fixed-Free 


Fig  2-21  Pressure  and  radial  velocity  vs  non-dimensional  radius  for  the  first  three 
free-free  modes  in  a  fluid  tube.  Solid  line  =  pressure;  dashed  line  =  velocity 
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Fig  2-22  Pressure  and  radial  velocity  vs  non-dimensional  radius  for  the  first  three 
free-fixed  modes  in  a  fluid  tube.  Solid  line  =  pressure,  dashed  line  =  velocity. 
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Fig  2-23  Pressure  and  radial  velocity  vs  non-dimensional  radius  for  the  first  three 
fixed-free  modes  in  a  fluid  tube.  Solid  line  =  pressure,  dashed  line  =  velocity. 
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Fig  2-24  Pressure  and  radial  velocity  vs  non-dimensional  radius  for  fixed-fixed 
modes  in  a  fluid  tube.  Plane  wave  mode  not  shown.  Solid  line  =  pressure, 
dashed  line  =  velocity. 
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Fig.  2-26  Normalized  specific  acoustic  impedance  of  external  fluid  as 
a  function  of  non-dimensional  frequency  ka  for  zero  axial  wavenumber  kz. 


(b)  Reactive  component  of  normalized  external  impedance 
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Figure  2-27  Normalilzed  specific  acoustic  impedance  of  external  fluid  as  a 
function  of  non-dimensional  frequency  for  four  values  of  axial  wavenumber. 
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(a)  Real  component  of  normalized  external  impedance 


(b)  Imaginary  component  of  normalized  external  impedance 


Fig.  3-1  Branch  1  (solid  line)  and  Branch  2  (dashed  line)  frequencies  for  thin  elastic 
shell  as  functions  of  axial  wavenumber.  Frequencies  in  kilohertz,  wavenumber  in 
1/meters 
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Fig.  3-2  Phase  velocity  of  Branch  1  (solid  line)  and  Branch  2  (dashed  line)  for  thin 
elastic  shell  as  functions  of  axial  wavenumber.  Velocities  in  meters/second, 
wavenumber  in  1 /meters 
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Fig.  3-3  Group  velocity  of  Branch  1  (solid  line)  and  Branch  2  (dashed  line)  for  thin 
elastic  shell  as  functions  of  axial  wavenumber.  Velocities  in  meters/second, 
wavenumber  in  1 /meters 


Fig.  3-4  Ratio  of  radial  to  axial  displacements  for  Branch  1  (solid  line)  and 
Branch  2  (dashed  line)  for  thin  elastic  shell  as  functions  of  axial  wavenumber. 
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1Q<  Fig.  3-5  (c)  RootQ3  Variation  with  Frequency 


Fig.  3-6  (a)  RootQI  Variation  in  Transition  Region 
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Magnitude  of  Complex  Roots 
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Fig  3-10  Magnitude  of  S  Roots  in  Transition  Zone 


Fig  3-11  (a)  Components  of  S3  in  Transition  Zone,  Beta 
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Fig  3-1 3  (a)  Magnitude  of  W(f,z)  vs  Distance  and  Frequency 
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Fig  3-14  (a)  Magnitude  of  U(f,z)  vs  Distance  and  Frequency 


Fig  3-14  (b)  20*!og  Magnitude  of  U(f,z)  vs  Distance  and  Frequency 


Fig  3-15  (a)  Magnitude  of  W(f,z)  vs  Distance  and  Frequency 


Fig  3-15  (b)  20*!og  Magnitude  of  W(f,z)  vs  Distance  and  Frequency 


Fig  3-16  (a)  Magnitude  of  U(f,z)  vs  Distance  and  Frequency 


Fig  3-16  (b)  Magnitude  of  U(f,z)  vs  Distance  and  Frequency 


Axial  Distance,  cm. 
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Fig  3-21  Magnitude  of  W(f,z)  vs  Distance  at  Fixed  Frequencies 


Fig  3-22  Magnitude  of  U(f,z)  vs  Distance  at  Fixed  Frequencies 
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Fig.  3-27  Spectogram  countour  plot  of  W(t)  at  Z  =  6  cm 
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Fig.  3-29  3-D  Spectogram  Contour  Plot  of  W(t)  at  Z  =  6  cm 
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Fig.  3-30  3-D  Spectogram  Contour  Plot  of  U(t)  at  Z  =  6  cm 
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Figure  4-2  Normalized  Spectral  Impedance  of  a  Fluid  Tube  in  an  External  Fluid  Impedance  vs 
Wavenumber  for  fluid  impedance  ratios  Rz  =  .5  (solid  curve),  1  (dotted  curve)  and  2  (dashed 
curve).  Sound  speed  ratio  Rc  =  2;  ratio  of  radii:  s  =  2 
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(d)  Imaginary  Component  of  Impedance 
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Fig.  5-2  Admittance  of  Disk  Near  First  Radial  Mode 


Fig.  5-3  Admittance  of  Disk  vs  Frequency 
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Fig.  5  —  6C b 3  Equivalent  Circuit  for  Measurements  in  Air 


Fig.  5-6CcD  Equivalent  Circuit  for  Dummy  Load 
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Fig.  5-8  Difference  between  computed  and  measured  modal  frequencies  (in  kHz)  vs  mode 
number  for  the  first  eight  modes.  Solid  curve  from  the  membrane  model,  dashed  curve 
represents  a  free-free  bar.  Calculated  values  based  on  wave  speed  of  3939  meters/second. 


Fig.  5-9  Standard  deviation  of  the  difference  between  measured  and  computed  modal 
frequencies  (in  Hz)  as  a  function  of  assumed  wave  speed.  Solid  curve  is  model,  dashed 
curve  is  free-free  bar 
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Fig.  5-10  Difference  between  computed  and  measured  modal  frequencies  vs  mode  number 
expressed  as  a  percentage  of  the  measured  frequency.  Model  results  indicated  by  V\  free-free 
bar  by 

Caclulated  modal  frequencies  based  on  extensional  wave  speed  of  3939  meters/second 
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Fig.  5-1 3  Spectogram  contour  plot  of  waveform  at  gage  C-3 


Fig.  5-14  Spectogram  contour  plot  of  waveform  at  gage  C-1 
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Fig.  5-1 7  3  D  Contour  plot  of  spectogram  of  waveform  at  gage  C-3 
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Fig.  5-1 8  3  D  Contour  plot  of  spectogram  of  waveform  at  gage  C-1 
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Fig.  5-1 9  3  D  Contour  plot  of  spectogram  of  waveform  at  gage  C-4 


Fig.  5-20  3  D  Contour  plot  of  spectogram  of  waveform  at  gage  C-2 
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Table  2-1 

Roots  of  the  characteristic  equation  for  free  wave  propagation  in  a  fluid  tube 
with  free-free  radial  boundaries,  for  four  values  of  radii  ratio  s.  Columns  labeled 
"Root( )"  show  the  exact  value  of  the  root  for  the  mode  indicated.  Columns 
labeled  "E  ( )"  show  the  error  between  the  asymptotic  estimate  and  the  exact 
root,  expressed  as  a  percent  of  the  exact  root. 


Mode 

Number 


s  =  2 


Mode  s  -  3  s  =  5 

Number 


m  Root3  E3  Root4  E4 


Table  2-2 


Roots  of  the  characteristic  equation  for  free  wave  propagation  in  a  fluid  tube 
with  free  boundary  at  inner  radius  and  fixed  boundary  at  outer  radius,  for  four 
values  of  radii  ratio  s. 


Mode 

Number 

S 

=  1.5 

s  =  2 

m 

Rootl 

m 

El 

m 

Root2 

m 

E2 

m 

i 

9.344 

0.865 

4.645 

1.455 

2 

15.66 

0.307 

7.814 

0.513 

3 

21.957 

0.156 

10.967 

0.26 

4 

28.248 

0.094 

14.115 

0.157 

5 

34.536 

0.063 

17.261 

0.105 

6 

40.822 

0.045 

20.405 

0.075 

7 

47.108 

0.034 

23.549 

0.056 

8 

53.393 

0.026 

26.692 

0.044 

9 

59.678 

0.021 

29.835 

0.035 

10 

65.962 

0.017 

32.977 

0.029 

Mode 

Number  s  =  3 

m  Root3  E3 


m  _ m 


1 

2.303 

2.326 

2 

3.895 

0.819 

3 

5.475 

0.416 

4 

7.051 

0.251 

5 

8.625 

0.168 

6 

10.198 

0.12 

7 

11.77 

0.09 

8 

13.342 

0.07 

9 

14.914 

0.056 

10 

16.486 

0.046 

s  =  5 


Root4  E4 

_ m  _ m 


1.139 

3.408 

1.939 

1.275 

2.731 

0.654 

3.52 

0.397 

4.308 

0.266 

5.095 

0.191 

5.882 

0.144 

6.668 

0.112 

7.455 

0.09 

8.241 

0.073 

Table  2-3 


Roots  of  the  characteristic  equation  for  free  wave  propagation  in  a  fluid  tube 
with  fixed  boundary  at  inner  radius  and  free  boundary  at  outer  radius,  for 
four  values  of  radii  ratio  s. 


Mode 

Number 


s=  1.5 


s  =  2 


m 

Rootl 

m 

El 

m 

Root2 

m 

E2 

m 

i 

3.401 

-7.636 

1.794 

- 12.439 

2 

9.52 

-0.998 

4.801 

-1.838 

15.766 

-0.367 

7.909 

-0.692 

4 

22.033 

-0.189 

11.035 

-0.357 

5 

28.307 

-0.114 

14.168 

-0.217 

6 

34.584 

-0.077 

17.304 

-0.146 

7 

40.863 

-0.055 

20.442 

-0.105 

8 

47.143 

-0.041 

23.58 

-0.079 

9 

53.424 

-0.032 

26.72 

-0.061 

10 

59.706 

-0.026 

29.86 

-0.049 

Mode 

Number 


m 

RooG 

E3 

Root4 

E4 

, . 

m 

m 

m 

m 

1 

0.959 

- 18.129 

0.515 

-23.706 

2 

2.435 

-3.238 

1.246 

-5.483 

3 

3.978 

-1.275 

2.009 

-2.264 

4 

5.535 

-0.669 

2.783 

- 1.227 

5 

7.098 

-0.409 

3.561 

-0.763 

6 

8.663 

-0.276 

4.342 

-0.519 

7 

10.23 

-0.198 

5.124 

-0.375 

8 

11.799 

-0.149 

5.907 

-  0.283 

9 

13.367 

-0.116 

6.691 

-0.221 

10 

14.936 

-0.093 

7.475 

-0.178 

Table  2-4 


Roots  of  the  characteristic  equation  for  free  wave  propagation  in  a  fluid  tube 
with  fixed  boundaries  at  both  inner  and  outer  radii,  for  four  values  of 
radii  ratio  s. 


Mode 

Number 

m 

T 

2 

3 

4 

5 

6 

7 

8 

9 

10 


s 

=  1.5 

S  =  2 

Rootl 

m 

El 

m 

Root2 

m 

E2 

m 

6.322 

-0.609 

3.196 

- 1.693 

12.586 

-0.157 

6.312 

-0.46 

18.863 

-0.07 

9.444 

-0.208 

25.143 

-0.04 

12.581 

-0.118 

31.424 

-0.025 

15.72 

-0.076 

37.706 

-0.018 

18.859 

-0.053 

43.988 

-0.013 

22 

-0.039 

50.27 

-0.01 

25.14 

-0.03 

56.553 

-0.008 

28.281 

-0.023 

62.836 

-0.006 

31.422 

-0.019 

Mode 

Number  s  =  3  s  =  5 


m 

Root3 

m 

E3 

m 

Root4 

m 

E4 

m 

i 

1.636 

-3.957 

0.847 

-7.253 

2 

3.178 

-1.161 

1.61 

-2.455 

3 

4.738 

-0.54 

2.385 

-1.211 

4 

6.303 

-0.309 

3.164 

-0.712 

5 

7.87 

-0.2 

3.945 

-0.466 

6 

9.438 

-0.139 

4.728 

-0.328 

7 

11.007 

-0.103 

5.511 

-0.243 

8 

12.576 

-0.079 

6.295 

-0.187 

9 

14.146 

-0.062 

7.079 

-0.148 

10 

15.716 

-0.051 

7.863 

-0.12 

Table  5-1  Published  Data  on  Brass 
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<u 

£ 

cc 

Lambda 

£ 

c 

2 

- 

9.45E+10 

1  08E+11 

1.08E+1 1 

0l+Ht9  8 

5.07E+10 

8.96E+10  1 

1.15E+1 1 

9.31E+10 

1.99E+10 

21.34 

1.19E+11 

Bulk 

Modulus 

£ 

a 

2 

- 

11+3611 

1.30E+1 1 

1.36E+11 

1.10E+11 

7.10E+10 

1.15E+11 

1.95E+1 1 

1.25E+11 

3.45E+10 

27.60 

2.02E+11 

C5 

<L> 

xi 
c n 

Modulus 

s 

a 

CO 

! 

3.67E+10 

3.80E+10 

3.80E+10 

3.53E+10 

3.04E+10 

3.78E+10 

4.03E+10 

01+399T 

2.90E+09 

7.91 

4.17E+10 

u* 

a 

o 

x; 

C/3 

Wavspd 

Cs,  m/s 

2080 

1948 

1948 

2030 

2050 

2120 

2173 

2050 

00 

t> 

3.78 

2211 1 

Thin  Bar 

Wavspd 

Cb,  m/s 

3430 

3492 

3498 

<N 

Tt* 

CO 

CO 

3305 

3485 

3597 

o 

CO 

CO 

o\ 

2.69 

3660 

C> 

03 

s 

Wavspd 

Cp,  m/s 

1 

r- 

r- 

'O 

co 

3759 

3765 

3574 

3463 

CO 

<N 

r- 

CO 

3872 

3690 

cn 

cs 

3.39 

3939 

Long. 

Wavspd 

1 

4450 

4644 

4700 

4280 

3830 

4430 

00 

r- 

4445 

299 

<N 

\d 

4868 

[  Poisson's 

Ratio 

VsO 

CO 

© 

0.37 

0.37 

0.35 

O 

CO 

o 

0.35 

CO 

O 

j 

0.35 

0.02 

6.60 

0.37 1 

So 

C 

D 

D 

Kg/cu  m 

8.50E+03 

8.53E+03 

CO 

O 

+ 

W 

o 

00 

8.56E+03 

CO 

O 

+ 

W 

o 

00 

8.42E+03 

8.53E+03 

8.45E+03 

1.48E+02 

1.75 

8.52E+03  | 

Young's 

Modulus 

|  N/sq  m 

1.00E+11 

1.04E+11 

1.04E+11 

9.56E+10 

8.85E+10 

1.02E+11 

_ 1.10E+1 1 

1.01E+1 1 

6.46E+09 

6.42 

1.14E+11 1 

Data  Source 

i 

Fany 

(Elmore  &  Heald 

Kinsler  &  Frey 

(Metrotek,  Brass 

Metro.  Half  Hard 

Metro,  Naval  Brass 

r 

1 

i 

&\ 

ml 

<L> 

c 

t: 

C3 

U_l 

i1 

ii 

_L 

AVERAGE 

t 

1 

j! 

2 

O 

H 

< 

> 

U 

2) 

P 

s> 

i 

< 

c 

2 

2 

Modal  Tests  | 

Table  5-2 


Comparison  of  measured  and  calculated  modal  frequencies  (in  Hz) 


Measured  modal 
frequencies 


Modal  frequencies 
calculated  from 
membrane  model 


Modal  frequencies 
calculated  for  free- 
free  bar 


meas 


'  506.38 

'506.834 

'  506.85 

1014.937 

1013.569 

1013.7 

1523.59 

1520.106 

1520.55 

2029.92 

2026.346 

2027.399 

2535.514 

*  model 

2532.185 

^nd 

2534.249 

i 

3042.169 

3037.519 

| 

3041.099 

3543.690 

3542.236 

3547.949 

4044.796 

.  4046.225 

:  4054.799 

The  calculated  modal  frequencies  were  computed  using  an  extensional 
wave  speed  of  3939  meters/second. 
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